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1  Summary 

The  drag  due  to  the  fuselage  and  hub  constitutes  a  major  portion  of  the  overall  drag  of  a 
rotorcraft.  Indeed,  “hub  drag"  (drag  from  the  rotor  hub,  swashplate,  blade  actuators,  hub 
pylons,  engine  nacelles  and  blade  shanks),  can  constitute  over  ten  percent  of  the  total  power 
required  on  typical  helicopter  designs,  and  over  forty  percent  of  parasite  drag.  On  modern 
compound,  high-speed  or  long-endurance  rotary-wing  vehicle  concepts,  increased  complexity 
may  result  in  the  combination  of  the  bluff  body  drag  contributors,  i.e.,  the  hub  and  fuselage, 
becoming  a  primary  limiter  of  performance.  Prior  to  the  commencement  of  this  effort,  the 
knowledge  of  where  the  state  of  the  art  in  predicting  this  large  drag  contribution  using 
computational  fluid  dynamics  (CFD)  methods  was  inadequate  to  predict  whether  a  design 
for  an  advanced  rotorcraft  will  be  close  to  actuality,  even  when  it  is  compared  to  one  or  two 
small-scale  model  experiments.  An  understanding  of  the  various  contributors  to  hub  drag 
and  their  scaling  with  Reynolds  number  and  other  flight  parameters,  is  critical  to  improving 
current  designs,  and  to  assessing  new  concepts. 

The  effort  specifically  involves  fundamental  research  that  systematically  investigated  the 
scaling  of  complex  nonlinear  effects  such  as  bluff  body  separation  and  rotational/stationary 
interference  effects  on  a  generic  hub  configuration  that  included  shanks,  blocks,  pitch  links 
and  scissors.  The  effort  proceeded  in  two  parts:  an  experimental  investigation  and  a  con¬ 
current  computational  investigation.  For  all  configurations,  except  where  experimental  wake 
data  collection  locations  were  not  known,  all  computations  were  computed  a  priori  to  ob¬ 
taining  the  experimental  data.  For  the  cases  where  the  simulations  could  not  be  completed 
before  the  experiments,  no  manipulation  of  the  simulations  from  the  original  protocol  es¬ 
tablished  occurred,  ensuring  a  priori  and  independent  correlation  of  the  experiment  and 
computations. 

The  primary  objective  of  computational  portion  of  this  project  was  to  reduce  the  uncer¬ 
tainty  in  current  numerical  predictions,  through  computational  investigations  that  leveraged 
a  Vertical  Lift  Consortium  (VLC)-funded  hub  drag  scaling  research  effort.  To  confirm  this 
objective,  correlations  are  performed  with  the  results  obtained  during  the  experimental  por¬ 
tion  of  this  project.  Experiments  were  performed  on  a  1/3.5  model  scale  experimental  model, 
which  included  data  from  load  cells,  as  well  as  selected  wake  measurements  using  a  combi¬ 
nation  of  particle  image  vclocimctry  (PIV)  and  hot-wire  anemometry. 

To  achieve  this  objective,  advanced  turbulence  methods  based  on  a  hybrid  unsteady 
Reynolds- averaged  Navier-Stokes/Large  Eddy  Simulations  (URANS/LES)  approach  was  ap¬ 
plied.  A  new  feature-based  anisotropic  grid  adaptation,  capable  of  grid  modification  on 
overset  meshes  has  been  developed,  validated  and  demonstrated  in  this  approach. 

The  objectives  of  the  computational  portion  of  this  research  were  all  met.  CFD  methods, 
using  the  approach  developed,  can  be  applied  on  a  single  baseline  grid  for  a  range  of  Reynolds 
numbers  and  operating  conditions.  This  eliminates  the  need  to  generate  grids  to  maintain 
accuracy  in  the  near  and  far  wake.  Excellent  correlation  with  experimental  performance  and 
wake  characteristics  were  obtained  using  a  priori  simulations  for  the  recommended  approach. 
Computational  analysis  was  also  used  to  identify  potential  errors  in  the  experimental  data. 

The  hub  characteristics  were  identified  with  computational  methods  for  both  full  assembly 
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and  deconstructions  of  the  hub.  From  these,  identification  of  where  theoretical  approxima¬ 
tions  can  be  applied  were  identified.  Interference  effects  were  also  quantified,  and  the  sources 
identified,  including  the  influence  of  a  fuselage.  Strouhal  shedding  and  wake  interactions  were 
computed  from  the  rich  data  provided  from  the  computational  simulations.  Important  ef¬ 
fects  of  scaling  for  both  static  and  rotating  hubs  were  quantified;  it  is  clear  that  there  are 
significant  differences  in  the  behavior  of  the  hub  and  the  component  interference  between 
model  and  full  helicopter  scales.  This  implies  the  critical  need  of  advanced  computational 
techniques  such  as  the  one  developed  here  wrhen  designing  and  analyzing  hubs.  Long-age 
wake  data  were  captured  and  correlated  with  experimental  data  from  other  sources  to  show 
that  prediction  and  analysis  of  empennage-based  aeroelastic  and  unsteady  aerodynamic  phe¬ 
nomena,  such  as  tail  buffet  and  “wag”  arc  within  reach  of  current  computational  resources. 
Turbulent  spectra  confirmed  the  ability  to  maintain  complex  wake  behavior  over  these  long 
periods. 

Computational  methods  have  now  reached  a  maturity  that,  if  care  is  taken  in  the  grid 
generation  and  turbulence  modeling,  it  may  be  more  cost-effective  and  accurate  to  design 
and  analyze  complex  hubs  directly  with  computational  methods  rather  than  reliance  on 
model  scale  experiments.  Advanced  turbulence  models  that  use  detached  or  large  eddy 
simulations  and  grid  adaptation  are  recommended  to  improve  the  performance  quantities 
and  the  unsteady  wake  characteristics. 
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2 


Nomenclature 


a 

b 


c 


C 

Cj1 

cv,ct 

U 

E 

£ 

T 

Fe,(~) 

Etoi 


h 

h° 

hi 

H 

H 


i,j  k 
I 
k 
k 

^ sgs 

L 

L 


M 


Mfip 

n 


n\,  n2 
N 


P 

Pr 


Qi 

Qcrit 

R 


Rec 


s 

S 

t 

T 


Speed  of  sound,  m/s 
Half  chord  length,  m 
Rotor  blade  chord,  m 
Pressure  coefficient 
Coarsening  factor 

Thrust  coefficient  (p^nfl)2-) 

Coefficients  used  in  LES  k  equations 
Time-averaged  fuselage  centerline  pressure  integral 
Total  energy,  ( E  =  CVT  +  ttfcita/2) 

RANS  statistical  filtering  operator 
LES  statistical  filtering  operator 

Adaptation  formulation  across  an  edge  based  on  a  solution  feature 

User  specified  adaptation  tolerance 

Specific  enthalpy 

Original  mesh  spacing,  grid  units 

New  mesh  spacing  request,  grid  units 

Hessian  matrix 

Additive  hybrid  blending  operator 
Cartesian  unit  vectors 
Adaptation  intensity 

Turbulent  kinetic  energy  or  reduced  frequency 
Subgrid  scale  turbulence  kinetic  energy 
Edge  length,  grid  units 
Leonards’  stress 

Anisotropic  adaptation  metric  or 
Stress  term  used  in  Cv 
Tip  Maeh  number 

Time  step  index  within  adaptation  window 

Node  (1  or  2)  of  given  edge 

Number  of  time  steps  within  adaptation  window 

Static  pressure,  N/m2 

Prandtl  number 

Heat-flux  vector  (<&  =  —ndT/dXi) 

Q-eriterion,  s-2 

Rotor  radius,  m,  or  gas  constant 

Reynolds  number  based  on  blade  chord  and  tip  speed  ^ 

Instantaneous  strain-rate  tensor 
Mean  strain-rate  tensor,  s_1 
Time,  s 

Rotor  thrust,  N 
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u,v, 

w 

= 

Instantaneous  local  velocity  components  in  (x,y,z),  m/s 

Uoo 

= 

Free  stream  velocity,  m/s 

vh 

= 

Effective  hot-wire  measured  velocity,  m/s 

Vtip 

= 

Rotor  tip  speed,  m/s 

w 

= 

Adaptation  time  window 

x,y, 

z 

— 

Streamwise,  spanwise  and  vertical  directions,  m 

X 

= 

Spatial  vector,  m 

y+ 

= 

Nondimensional  normal  distance  of  the  first  cell  from  the  wall 

Greek  Symbols 

Pis, 

Pic 

= 

Lateral  and  longitudinal  first  flapping  harmonics 

= 

Kronecker  delta 

A 

= 

Characteristic  length  for  grid  scale 

A 

= 

Characteristic  length  for  test  scale 

7 

= 

Ratio  of  specific  heats 

K 

— 

Coefficient  of  heat  conductivity 

P 

= 

Rotor  advance  ratio  or  coefficient  of  molecular  viscosity,  N-s/m 

hr 

= 

Eddy  viscosity  nondimensionalized  by  // 

u 

= 

Kinematic  viscosity,  m2/s 

P 

= 

Density,  kg/m3 

fp 

= 

Blade  azimuth  angle,  0 

T 

= 

Viscous  stress  tensor, 

n 

= 

Rotation  rate  tensor,  s_1 

M 

— 

Vorticity  magnitude,  s_1 

u 

= 

Specific  turbulence  dissipation  rate  or  vorticity  magnitude,  s_1 

n 

— 

Rotational  velocity  of  rotor,  rad/s 

Mathematics  Symbols 
d  =  Partial  derivative  operator 

(•)  =  Averaged  value 

(•)  =  Favrc-avcragcd  value 

(•)  —  Variable  filtered  at  test  level 

(.)SGS  =  Subgrid  scale 

(.)'  =  Fluctuating  value 

(.)"  =  Fluctuating  Favre-averaged  variable  value 

(.)t  —  Turbulent  variable 

( •  )oo  =  Freestream  value 
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3 


Background 


3.1  Problem  Motivation 

Parasite  drag  on  rotorcraft  can  become  a  crucial  factor  in  forward  flight  especially  during 
high  speed  flight1,2  and  can  limit  the  range,  maximum  speed,  and  payload  of  the  vehicle. 
Reduction  in  parasite  drag  can  improve  vehicle  stability  and  control1  and  significantly  de¬ 
crease  vibrational  and  blade  loads  to  reduce  vehicle  weight  and  extend  the  rotor  blade  life.4 
Hub  assemblies  for  single  main  rotor  helicopters  can  contribute  nearly  25%-30%  of  vehicle 
parasite  drag,  while  hub  assemblies  for  coaxial  rotors,  such  as  those  on  the  XH-59,5  can 
contribute  as  much  as  50%  of  the  parasite  drag.  The  complex  and  highly  turbulent  rotor 
and  hub  wake  is  convcctcd  to  the  empennage  and  results  in  acroclastic  behavior  such  as  tail 
buffet  and  aerodynamic  behavior  knowm  as  tail  “wag” .  These  phenomena  result  in  fatigue 
and  reduced  handling  qualities.  Therefore,  reduction  of  the  drag  in  the  design  of  hub  systems 
can  be  critical  to  the  success  of  high-speed  rotorcraft  design.  In  order  to  achieve  these  goals, 
the  drag  sources  associated  with  complex  hub  designs  must  be  thoroughly  investigated,  well 
understood,  and  reliably  predicted  using  analysis  tools. 

There  is  a  large  literature  base  on  flight  to  wind  tunnel  correlation  of  drag  for  fixed  wing 
vehicles6,7  preceding  the  establishment  of  cryogenic  wind  tunnels  capable  of  separately  or  si¬ 
multaneously  simulating  flight  Reynolds  and  Mach  numbers  with  sub-scale  models.  However, 
rotary  wing  designers  must  deal  with  acrodynamicallv  dirty  appendages  and  interference  of 
rotating  components,  where  the  issue  of  drag  is  much  more  complex  and  significant. 

Industry  estimates  that  a  thirty  percent  reduction  in  hub  drag  below7  the  current  empirical 
gross  weight  correlation  is  necessary  for  next  generation  vehicles.  Present-day  computational 
fluid  dynamics  (CFD)  codes,  especially  commercial  codes  suitable  for  vehicle  flow  field-based 
drag  analysis,  are  quite  inadequate  for  drag  prediction  even  at  the  conceptual  design  stage. 
This  point  is  reinforced  by  a  recent  NASA/DFLR  blind-prediction  exercise  on  a  relatively 
clean  cargo  aircraft  configuration  wThere  some  wing-root  separation  occurred  at  a  lift  coeffi¬ 
cient  of  0.5.  This  showed  a  wdde  scatter  (over  23%  uncertainty)  in  the  drag  polar  predictions 
from  34  different  CFD  researchers,8,9  wThcn  compared  wTit,h  experimental  results.  The  cited 
causes  for  uncertainties  are  the  lack  of  accurate  methods  for  the  prediction  of  separation  on 
smooth  aerodynamic  surfaces,  as  well  as  the  application  of  RANS-based  turbulence  mod¬ 
els  that  are  pre-tuned  to  a  select  set  of  test  cases  and  may  not  suitable  for  the  physics  of 
the  configuration.  While  RANS  turbulence  models  typically  provide  acceptable  prediction 
of  separation  in  the  case  of  airfoil  leading  edge  stall,10  trailing  edge  stall  mechanisms  con¬ 
tinue  to  elude  capture  via  statistical  modeling.11  In  the  case  of  rotorcraft  hub  drag,  several 
mechanisms  of  separation  may  be  present  leading  to  poor  CFD  results,12  and  the  problem 
is  complicated  by  a  number  of  superimposed  nonlinear  concerns  including,  but  not  limited 
to  incompressible  flight  regime,  rotational  dowmwash  from  the  rotor,  and  transition  from 
laminar  to  turbulent  flow7.  Different  CFD  methods  are  unable  to  consistently  capture  these 
flows  when  applying  the  same  RANS  turbulence  closures  and  inadequately  comprehend  the 
impact  of  sealing  from  model  to  full  scale  Reynolds  numbers. 

These  inconsistencies  and  gaps  in  the  known  physics  are  not  limited  to  numerical  simu- 
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lations.  Keys  and  Rosenstein11  showed  that  at  low  advance  ratio  (/i),  an  increase  in  average 
shank  dynamic  pressure  explains  the  drag  increase.  Wind  tunnel  tests  on  compound  rotor- 
craft  show  large  uncertainty  in  source  and  magnitude  of  “interference55  drag  (i.c.,  discrepancy 
between  the  sum  of  component  drag  values  and  the  integrated  powered  vehicle  drag),  which 
must  be  understood  to  apply  hub  drag  reductions  consistently,  especially  during  the  reliance 
on  smaller  model-scale  wind  tunnel  tests. 

A  typical  hub  is  a  plethora  of  interacting  bluff  bodies,  which  are  associated  with  poor 
performance  characteristics  and  higher  drag.  The  primary  drag  component  in  these  bluff 
bodies  is  due  to  flow  separation  (pressure  drag)  rather  than  viscous  effects  (friction  drag). 
Bluff  body  wakes  associated  affect  the  performance  of  both  commercial  and  military  air  vehi¬ 
cles,13  in  particular  impacting  tail  component  fatigue  and  handling  qualities.  Accurate  hub 
drag  predictions  are  inhibited  by  many  differing  complex  flow  interactions.  Some  hub  wake 
characteristics  include  periodic  forcing  and  vortex  interactions.  In  many  cases,  these  inter¬ 
actions  originate  from  fine  structures  such  as  tubes,  wires  and  linkages.  The  identification  of 
the  influence  of  the  full  assembly  and  component  Reynolds  numbers  can  not  be  minimized; 
over  the  Reynolds  number  range  of  interest,  distinct  bifurcations  produce  significant  and 
measurable  differences  in  the  wake  flow  field.14 

Much  of  the  experimental  research  of  hub  drag  is  directed  toward  improving  drag  char¬ 
acteristics  of  current  hub  designs  by  the  addition  of  fairings.  Fairing  designs  have  been 
explored  by  Sikorsky  and  others  to  reduce  flow  separation  and  interference  drag  between 
the  hub  and  fuselage. 1  To  date,  frontal  swept  area  of  the  hub  design  has  been  the  lead¬ 
ing  parameter  tied  to  hub  drag,15  therefore  the  fairing  of  an  existing  hub  design  does  not 
address  the  issue  directly.  This  is  especially  true  for  articulated  hubs,  where  empty  space 
is  required  for  the  control  hinges,  resulting  in  channel  flows  that  interact  with  cylindrical 
components.  While  empirically  corrected  analytic  estimates  have  been  developed  to  predict 
hub  drag  based  on  frontal  area,  there  is  no  consistent  trend  when  accounting  for  interference 
effect  and  frontal  swept  area.15  Considerations  for  hub  displacement  from  the  fuselage  have 
been  made,  weighing  the  effects  of  increased  frontal  area  to  decreased  interference.15 

Based  on  these  issues,  this  research  was  designed  to  aid  in  meeting  the  needs  for  both 
commercial  and  military  next  generation  rotoreraft  through  the  reduction  of  drag  through 
the  methodologies  and  results  developed  herein.  Drag  reduction  is  also  a  major  goal  of  the 
DoD  thrusts.  This  effort  also  is  applicable  to  future  vehicle  derivatives  and  joint  transport 
rotoreraft  hub  drag  reduction.  It  can  be  applied  directly  on  11-92,  VII-92,  S-92  Growth. 
S-76,  X2  and  V22  vehicles. 

3.2  Computational  State  of  the  Art  for  Rotating  Hub  Analysis 

The  portion  of  the  research  undertaken  in  this  ONR-funded  project  reported  here  is  on  the 
computational  prediction  of  complex  hub  performance  and  their  turbulent  wakes.  Recently, 
there  have  been  several  fundamental  experimental  studies  of  hub  drag,  accompanied  by  state- 
of-the-art  computational  fluid  dynamics  (CFD)  predictions  of  these  complex  flows.5, 17-21 
These  studies  have  all  focused  on  models  that  are  a  fraction  (1/5  -  1/4  scale)  of  the  full- 
scale  rotor  hub.  Thus,  the  scaling  of  these  complex  hub  systems,  including  rotational  and 
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interference  effects,  must  also  be  understood. 

Wake  et  al.5  investigated  1/4-scale  faired  hubs  for  the  X2  Technology™  Demonstrator 
aircraft  using  an  unstructured  computational  solver.  These  simpler  faired  elliptical  geome¬ 
tries  can  prove  to  be  challenging  because  of  the  difficulty  of  predicting  separation  (and  po¬ 
tentially  transition).  The  focus  of  their  effort  was  to  investigate  the  impact  of  aerodynamic 
fairings  on  drag  for  the  dual  hub  configuration.  Using  a  grid  refinement  study,  they  were 
able  to  obtain  agreement  with  experiment  within  15%  for  their  tetrahedral  grids,  within  3% 
for  their  hexahedral  grids,  and  they  matched  within  8%  the  drag  estimates  obtained  by  two 
experimental  studies  for  various  configuration  changes.  The  configuration  that  they  analyzed 
was  static  and  did  not  include  components  such  as  root  stubs  or  hardware  in  the  analysis. 

A  follow-on  study  in  2011  by  Sikorsky1,  using  another  unstructured  method  confirmed  the 
overall  findings  of  their  initial  undertaking.  Similar  issues  have  been  recently  encountered 
during  another  experimental-computational  collaboration.21 

Bridgeman  and  Lancaster1*  have  studied  a  1/5-scale  Bell  rotor  hub  and  fuselage  both  ex¬ 
perimentally  and  computationally.  Using  an  extensive  grid  independence  study,  they  found 
that  total  drag  predictions  within  5% -10%  of  experimental  values  using  an  unstructured 
state-of-the-art  solver  (FUN3D)  could  be  achieved  for  the  non-rotating  hub-fuselage  config¬ 
uration,  depending  on  the  grid  resolution.  Details  such  as  hardware,  pitch  links  and  root 
stubs  were  included  in  the  eomputational  model;  a  breakdown  of  the  individual  contribu¬ 
tions  of  these  components  was  not  part  of  the  focus  of  this  work.  A  follow-on  effort  citcbcll2 
indicated  that  comparable  results  could  be  achieved  by  any  of  the  unstructured  solvers  that 
were  evaluated,  including  the  solver  utilized  in  this  study.  An  extension  to  the  methodology 
for  rotating  hubs  was  performed  by  Hill  and  Louis.20 

Reich  et  al.21  studied  a  notional  rotor  hub  at  1/3  and  2/3  scale  Reynolds  number  both 
computationally  and  experimentally  at  a  single  advance  ratio  of  0.2  in  the  Penn  State  Uni¬ 
versity  water  tunnel.  Measurements  from  experiment  included  total  hub  drag  and  wake 
diagnostics.  Computations  with  a  commercial  code  were  also  performed.  They  found  that 
the  most  prominent  wake  structures  were  from  the  two-per-revolution  (scissors)  and  four- 
pcr-revolution  (main  hub  arms).  Both  of  these  vortical  structures  persisted  far  downstream 
of  the  hub,  with  the  four-per-rev  structure  dissipating  faster  than  the  two-per-rev  struc¬ 
ture.  Correlations  with  computations  were  not  as  accurate  as  the  analyses  performed  by 
Bridgeman  and  Lancaster,1*  10  and  the  computations  did  not  include  turbulent  wake  spectra 
correlations  with  experiment. 

For  most  of  these  simulations,  existing  unsteady  Reynolds- Averaged  Navier-Stokes  (URANS) 
turbulence  models  or  Detached  Eddy  Simulations  (DES)  were  utilized  to  obtain  the  numer¬ 
ical  simulation.  Studies  show  that  these  URANS  models  are  failing  when  the  flow  field 
includes  viscous-dommated  features,  such  as  separation  or  vortex  interactions.  This  poor 
correlation  is  not  surprising,  as  the  URANS  models  arc  statistical  approximations  of  the 
turbulence  scales.  Large  eddy  simulation  (LES)  to  direct  numerical  simulation  (DNS)  are 
needed  to  resolve  the  most  energetic  scales  that  dictate  the  behavior  of  vortex  shedding 
and  interaction  encountered  in  these  flows.  LES  is  capable  of  capturing  the  larger  eddies 
and  models  the  smallest  or  subgrid-scale  (sgs)  eddies,  permitting  coarser  grids  than  DNS, 
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which  captures  eddies  of  all  sizes.  Unfortunately,  classic  DNS  and  LES  require  significantly 
larger  grids  and  more  restrictive  temporal  step  size  to  capture  the  most  energetic  turbu¬ 
lent  features,  which  scale  nonlincarly  with  increasing  Reynolds  number.  Hybrid  turbulence 
techniques  have  become  an  area  of  interest  to  close  the  gap  between  RANS  and  LES/ DNS. 
Hybrid  techniques  in  legacy  CFD  codes  should  capture  the  largest  scales  wherever  grid  res¬ 
olution  is  sufficient  to  support  LES.  Thus  even  coarse  grids,  more  suitable  for  Very-Large 
LES  (VLES),  should  be  more  representative  of  the  flow  physics.  Of  these  hybrid  meth¬ 
ods,  the  most  well  known  is  Detached  Eddy  Simulation  (DES)  typiealfy  applied  with  the 
Mcntcr  k-w  SST  or  Spalart-Allmaras  URANS  turbulence  models.  LES  simulations  occur  in 
separated  flows,  following  the  Smagorinsky  eddy  viscosity  assumption,  while  RANS  turbu¬ 
lence  modeling  is  applied  in  attached  flows.  Two  recent  Georgia  Tech  hybrid  RANS/LES  and 
VLES/LES  turbulence  simulation  techniques  have  been  implemented  into  the  computational 
methodologies,  OVERFLOW  and  FUN3D  with  good  success  on  airfoils  and  rotors. 

The  approach  in  this  effort  is  to  integrate  a  priori  computational  predictions  with  exper¬ 
iments,  so  that  the  basic  flow  phenomena  are  carefully  understood  and  the  predictions  are 
validated  against  experiment  at  the  flow  field  phenomenon  level,  not  just  at  the  integrated 
performance  level.  This  provides  better  confidence  as  the  predictions  arc  extended  to  other 
configurations  and  conditions  than  those  in  the  validation  experiment. 

3.3  Computational  State  of  the  Art  for  Complex  Unsteady  Tur¬ 
bulent  Wakes 

The  resolution  of  unsteady  wake  features  is  also  essential  for  a  multitude  of  aeroelastic  appli¬ 
cations  pertinent  to  rotorcraft.  These  include,  but  arc  not  limited  to  modeling  of  aeroelastic 
rotor  blades  in  forward  flight,  rotor-fuselage  or  rotor-rotor  interaction,  helicopter-ship  in¬ 
teraction,  and  tail  buffet.  Other  applications  of  interest  can  include  wing-store  separation, 
ship-wake  interactions,  and  wind  turbines.  The  current  most  popular  computational  fluid 
dynamics  (CFD)  approach  to  resolve  these  multiple  reference  frame  applications  is  via  overset 
grids,  where  the  moving  body  or  component  meshes  are  generally  highly  refined  and  overset 
on  one  or  more  static  background  grids.22-26  Despite  the  ability  of  unstructured  overset 
methods  to  model  dynamic  bodies,  it  does  not  address  the  issue  of  numerical  dissipation 
that  can  result  in  inaccuracy  of  the  wake  physics.26,27 

Feature-based  grid  adaptation  for  unsteady  problems  has  been  applied  on  single  grids 
using  various  methodologies.  Accurate  predictions  of  hovering  rotors  in  a  single  rotating 
adaptive  mesh  have  been  performed  by  several  researchers.28-31  However,  these  scenarios 
can  not  be  immediately  applied  to  the  prediction  of  rotors  in  forward  flight,  where  adapta¬ 
tion  is  needed  in  both  the  background  inertial  reference  frame  and  the  near-bodv  rotating 
frame.  Further,  the  interaction  of  rotors  with  non-moving  bodies  such  as  fuselages,  wind- 
tunnel  struts,  and  other  configurational  components  also  require  moving-grid  capability  to 
simulate  multiple  motion  frames.  As  an  alternative  to  the  overset  configuration,  Park  and 
Kwon32  have  demonstrated  an  unstructured  sliding  mesh  approach  where  a  rotating  grid 
communicates  with  a  stationary  background  grid.  Here,  articulation  of  the  rotor  blades  was 
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made  possible  using  grid  deformation  based  on  a  spring  analogy,  and  also  applies  feature- 
based  grid  adaptation.  Another  non-overset  based  approach  has  been  described  by  Cavallo 
ct  al.,33  which  uses  unstructured  grid  movement  and  deformation  to  enable  moving  body 
adaptation. 

Past  research  efforts  in  overset  adaptation  have  in  many  instances  relied  on  an  off-body 
(background)  Cartesian  grid-based  adaptive  capability.  Mcakin34  presented  a  grid  compo¬ 
nent  grouping  algorithm  with  overset  structured  grids  using  a  method  of  adaptive  spatial 
partitioning  and  refinement  and  applied  it  to  background  Cartesian  grids.  Variations  of  this 
technique  have  been  subsequently  demonstrated  by  Henshaw  and  Schwcndcman35  and  Kan- 
nan  and  Wang.36  Canonne  et  al.3'  used  an  overset  structured  cylindrical  grid  topology  to 
simulate  rotor  motion  in  hover  where  the  background  grid  is  adapted. 

Hybrid-solver  developmental  efforts  have  focused  on  rotor  methodologies  where  two  sepa¬ 
rate  solvers  are  applied  in  the  near  body  and  background  regions,  respectively.  Duque  et  al.38 
have  employed  a  structured  near-body  and  unstructured  wake  grid  approach  to  evaluate  ro¬ 
tors.  Here,  isotropic  adaptation  was  applied  on  the  unstructured  background  grid  operating 
in  a  non-inertial  reference  frame.  This  work  recommended  the  use  of  anisotropic  adaptation 
to  accurately  capture  inherently  anisotropic  phenomena  such  as  tip-vortices  and  to  exploit 
computationally  efficiency  of  this  grid  adaptation  technique.  Park  and  Darmofal30  intro¬ 
duced  a  parallel  anisotropic  adaptation  capability  for  non-overset  tetrahedral  grids.  This 
technique  has  been  applied  to  investigate  several  applications  such  as  sonic-boom  propaga¬ 
tion,40  41  viscous  transonic  drag  prediction,40  and  re-entry  vehicle  configurations.42  Recent 
development  has  focused  on  the  coupling  of  a  body-fitted  unstructured  solver  with  a  high- 
ordcr  Cartesian  solver  to  propagate  the  wake  in  the  mid  and  far  fields.  Sankaran  ct  al.43 
and  Wissink  et  al. 44,45  have  successfully  implemented  automatic  mesh  refinement  (AMR)  in 
the  Cartesian  background  solver  to  resolve  the  wake  based  on  flow  field  features. 

The  effect  of  time-dependency  on  the  flow  field  is  an  essential  aspect  in  applying  grid 
adaptation  to  study  dynamic  moving  bodies  and  their  wake  structure.  Researchers33  37  have 
shown  that  adapting  the  solution  at  a  given  frequency  (based  on  flow  time)  have  proven 
effective,  with  increasing  frequency  yielding  higher  accuracy.  Investigations  involving  off- 
body  Cartesian-based  adaptive  mesh  refinement35, 36,43-45  have  extended  this  rationale  to 
adapt  the  solution  a  frequency  comparable  to  that  of  the  solver  time  step,  allowing  for  a 
coupled  adaptive  flow  solver.  However,  this  capability  is  not  computationally  efficient  for 
tetrahedral  unstructured-based  methods  as  they  do  not  have  the  advantage  of  octree  data 
structures  to  provide  for  the  adaptive  mechanics. 

An  alternative  approach  is  thus  required  which  addresses  the  time-dependency  issue 
without  the  computational  overhead  of  frequent  adaptation.  Kang  and  Kwon29  present  an 
adaptation  technique  that  detects  local  maxima  of  a  vortex  core  every  five  degrees  and  using 
a  3-D  parabolic  blended  curve  to  represent  the  vortex  core  path.  Park  and  Kwon32  describe 
a  '  quasi-unsteady  ’  adaptive  procedure  for  rotors  in  forward  flight  based  on  a  time  period  or 
window  dependent  on  the  blade  passing  frequency.  Cells  satisfying  an  adaptation  indicator 
arc  marked  at  each  time  step  within  the  window  and  adaptation  is  performed  for  those  cells  at 
the  conclusion  of  each  window.  Sterenborg  et  al.46  describes  a  similar  technique  and  applied 
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it  for  a  fluid-structure  interaction  investigation.  The  extension  of  a  time-dependent  feature- 
based  adaptation  methodology  for  anisotropic  grids  involving  dynamic  bodies  is  delineated 
by  Alauzct  and  Oliver.47 

These  successful  implementations  of  grid  adaptation  provide  impetus  for  further  investi¬ 
gation  of  this  approach  for  rotorcraft  applications  involving  multiple  grids.  What  is  common 
across  these  prior  overset-based  efforts  is  that  the  adaptation  is  restricted  to  the  off-body 
background  meshes.  Since  vorticity  originates  on  a  viscous  surface  where  the  near-bodv  grid 
is  employed,  the  full  capability  of  the  adaptation  cannot  be  exploited  unless  the  adaptation 
can  occur  across  the  meshes.  In  addition,  for  rotorcraft  configurations,  a  timc-dcpcndcnt 
strategy  is  necessary  to  permit  accurate  and  efficient  simulation  of  the  wake  growth. 
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4  Computational  Objectives  and  Goals 

The  long  range  computational  goals  of  this  project  is  demonstrate  the  capability  of  com¬ 
putational  methods  (CFD)  for  use  on  hub  and  hub/fuselage  design  and  analysis  with  high 
confidence.  In  addition,  the  ability  to  predict  the  complex  wake  structure  for  near  and  far 
wakes  is  needed  to  address  other  design  and  analysis  issues  on  the  helicopter. 

4.1  Objectives 

The  objectives  of  the  computational  portion  of  the  project,  as  outlined  in  the  original  pro¬ 
posal  are: 

1.  Isolate  and  quantify  the  different  sources  of  hub/pylon/naccllc  drag. 

2.  Tighten  the  tolerances  of  upper-bound  empirically-based  computations  suitable  for 
conceptual  design,  by  reference  to  the  basic  experiments. 

3.  Refine  the  capability  to  predict  drag  from  first  principles  through  computational  aero¬ 
dynamics,  suitable  for  the  preliminary  design  stage  and  beyond. 

4.  In  the  process,  advance  the  state  of  knowledge  on  predicting  flows  around  complex 
configurations  involving  flowr  separation  from  fixed  and  smooth  surfaces,  as  well  as  a 
wide  range  of  Rcvnolds  numbers. 
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5  Computational  Methodology  and  Development 

5.1  Baseline  Computational  Methodology 

The  baseline  grid  and  solver  methodologies,  prior  to  modification,  are  described  here. 

5.1.1  Grid  Generation  Approach 

NASA's  VGRIDns  grid  generation  software48  was  used  to  generate  either  fully  tetrahedral  or 
mixed  clement  meshes  for  configurations  studied.  To  perform  a  rotating  simulation,  meshes 
were  overset  to  form  a  composite  grid  using  SUGGAR++,49  the  overset  assembly  code 
documented  in  Section  5.1.2.  The  unstructured  overset  grid  methodology  allows  for  mesh 
refinement  about  complex  shapes,  and  limits  the  number  of  overset  meshes  needed  to  perform 
the  computation.  For  most  of  the  configurations  discussed  in  this  report,  comparisons  were 
made  with  wind  tunnel  experiments,  requiring  the  background  grid  to  represent  the  tunnel 
test  section. 

Where  ever  possible,  the  viscous  boundary  layers  of  the  components  were  modeled  using 
prismatic  elements,  while  the  mesh  surrounding  the  configuration  was  comprised  of  tetrahe¬ 
dral  elements.  This  mixture  of  elements  has  been  shown  to  accurately  capture  the  behavior 
of  the  viscous  boundary  layer  while  minimizing  the  cells  needed  in  the  external  flow  field.50 
Based  on  prior  experience  with  accurately  predicting  viscous  effects,51  between  30-50  nor¬ 
mal  mesh  layers  were  used  to  model  the  boundary  layers.  The  baseline  mesh  external  to  the 
boundary  layer  was  generated  with  general  wake  refinement,  such  that  a  new  grid  adaptation 
strategy  (Section  5.2)  was  used  to  modify  the  grid  based  on  the  specific  fiowfield. 

Details  of  the  hub  model  (including  the  bolt  nuts,  holes),  the  mount  and  the  wind  tunnel 
test  section  were  included  to  correlate  with  the  wind  tunnel  experiments.52  The  grid  fidelity 
are  further  described  in  each  individual  geometry  section. 

5.1.2  Numerical  Solver  Description 

FUN3D,  NASA’s  unstructured  RANS  solver,  was  selected  as  the  flow  solver  with  which  to 
demonstrate  the  new  adaptation  strategy.  FUN3D  utilizes  an  implicit,  node-based  finite 
volume  scheme  to  resolve  the  RANS  equations  on  unstructured,  mixed-topology  grids.53 
Both  compressible  and  incompressible54  Mach  regime  capabilities  arc  available  in  the  flow 
solver.  Time-accuracy  is  achieved  using  a  second-order  backward  differentiation  formula 
(BDF).  Roe’s  flux  difference  splitting  scheme55  is  used  compute  the  inviscid  fluxes,  while 
an  equivalent  central  difference  approximation  is  utilized  to  resolve  the  viscous  fluxes.  A 
Gauss-Seidel  strategy  is  used  to  solve  the  resulting  linear  system  of  equations.  FUN3D  has 
available  a  plethora  of  turbulence  methods,  of  which  Menter’s  ku-SST56  two-equation  model 
and  a  single-equation  turbulent  kinetic  energy  large  eddy  simulation  (LES)  model57,58  were 
applied  in  this  effort. 

Overset  functionality,  which  is  essential  to  simulate  grid  motion  of  a  rotating  hub  in 
the  presence  of  a  non-rotating  test  facility  or  fuselage,  is  achieved  via  two  libraries,  SUG- 
GAR++  and  DiRTlib.  The  Structured,  Unstructured,  and  Generalized  overset  Grid  As- 
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Figure  1:  Mean  pressure  coefficient  along  the  centerline  for  the  three-dimensional  circular 
cylinder  on  several  grids  with  varying  turbulence  methods.  LES  data  is  from  Ref.  65,  and 
experimental  data  is  from  Ref.  66.  From  63. 


sembleR  (SUGGAR++)  code  provides  for  Chimera  of  overset  capability  for  structured,  un¬ 
structured,  and  Cartesian  grid  topologies.49,59  Its  application  is  primarily  foi  composite  grid 
assembly  for  moving  body  simulations.  The  code  is  designed  to  perform  as  either  a  library  or 
as  a  stand-alone  executable  catering  to  various  grid  formats.  Composite  grids  can  be  created 
for  both  node  and  cell  centered  flow  solvers.  The  code  creates  the  composite  grid  along  with 
a  domain  connectivity  information  (DCI)  file.  The  I)CI  file  contains  information  of  the  as¬ 
sembly  of  multiple  overlapping  grids  via  interpolation.  It  lists  the  locations  of  fringe  points, 
where  interpolation  is  required,  as  well  as  donor  points  that  provide  interpolation  weights 
and  the  interpolated  values  at  the  fringe  locations.  Additionally,  hole-cutting  information 
and  orphans  (points  lacking  interpolation  data)  are  documented  in  this  file. 

FUN3D  s  overset  capabilities  have  been  successfully  applied  to  compressible  and  incom¬ 
pressible  rotorcraft  applications. 60~b2  In  such  simulations,  the  background  grid,  which  con¬ 
sists  of  the  fuselage  and  other  wind-tunnel  static  geometries  up  to  the  far-held  boundaries, 
is  assembled  with  finer  ncar-body  grids  for  each  of  the  moving  rotor  blades. 

5.1.3  Advanced  Turbulence  Modeling:  Hybrid  URANS-LES  Formulation 

A  hybrid  unsteady  Reynolds- Averaged  Navier-Stokes  and  large  eddy  simulation  (URANS/LES) 
turbulence  approach  (GT-HRLES)  within  FUN3D63  wTas  used  to  obtain  time-accurate  evalu¬ 
ations  of  the  unsteady  hub  wake.  This  turbulence  approach  has  been  showm  to  capture  both 
attached  and  separated  flow’s  more  accurately  than  URANS  alone.63,64  In  particular,  the 
separation  location,  wake  properties  and  performance  quantities  of  cylinders  is  comparable 
to  LES  results,  and  within  experimental  error  bounds.63  These  are  illustrated  in  Fig.  1  and 
Tabic  1. 
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Nz 

Turbulence 

model 

Mean 

cD 

Strouhal 

no. 

Separation 

location 

- 

- 

Experiment67 

0.99  ±0.05 

0.215  ±  0.00568 

86  ±2° 

101 

4  D 

kw-i SST 

1.456 

0.213 

98.4° 

2 

- 

DES 

1.5 

0.25 

86.8° 

51 

4  D 

HRLES 

0.971 

0.216 

85.8° 

48 

7 tD 

LES65 

1.04 

0.210 

88.0° 

Table  1:  Predicted  circular  cylinder  characteristics  for  various  turbulence  methods  and  grids 
with  number  of  spanwisc  planes,  Nz ,  and  spanwisc  extent,  Z.  Separation  location  is  given 
in  degrees  of  azimuth  from  the  leading  edge  stagnation  point.  From  Ref.  63. 


The  hybrid  URANS-LES  approach57,58  is  a  modeling  strategy  that  extends  a  URANS 
approach,  represented  by  the  statistical  operator  £  and  the  LES  approach  represented  on 
an  explicit  or  implicit  filtering  operator  T .  The  approach  relies  on  the  introduction  of  an 
additive  hybrid  operator  Ti  where  Tt  —  kpT  +  kp£  with  a  constraint  of  kp  +  kp  =  1  for 
two  blending  factors,  kp  and  kp ,  which  may  be  functions  of  space  and  time.  The  constraint 
kp  +  kp  =  1  is  very  important  as  it  imposes  the  preservation  of  the  constants,  which  is  a 
fundamental  property  necessary  for  an  averaging  operator.  In  this  effort,  the  URANS  model 
is  the  two-equation  Menter  Axj-SST  (shear  stress  transport)  model69  and  the  LES  approach 
solves  the  k  model.  In  the  current  hvbrid  RANS-LES  approach,  kp  =  1  —  kp ,  so  that 

H  =  (l-kE)T  +  kE£. 

Time  or  Reynolds-averaging  provides  the  typical  statistical  operation  for  the  URANS 
equations,  but  an  LES  approach  requires  a  mass-time  averaging  or  Favre  filtering.  The 
difference  between  these  two  averaging  or  filtering  approaches  is  key  to  resolution  of  the 
turbulence  in  each  approach.  In  URANS,  the  statistical  operator,  £,  resolves  the  turbulent 
velocity,  Ui(x,  t ),  into  fluctuating  and  mean  components,  t^(x,  t)  =  u'(x,  t)  +Uj(x,  £),  where, 
over  a  sufficient  time  period  (t  »  A'),  the  mean  of  the  fluctuating  velocity,  u'(x)  =  0. 
In  LES,  the  spatial  filtering  operator  T  separates  the  turbulent  velocity  into  a  large  scale 
component  and  a  small  or  sub  grid  scale  (SGS)  component.  The  SGS  component,  u',  differs 
from  the  Reynolds  decomposition  in  that  the  filtered  value  of  the  SGS  is  not  zero  ul  ^  0. 
Moreover,  the  LES  filtered  field  is  not  idempotent  (ii  /  ii),  as  is  the  assumption  for  the 
RANS  statistical  averaging.  The  fluctuating  LES  filtered  terms  arc  denoted  with  double 
prime  notation  (u"(x,  t))  to  differentiate  the  two  quantities. 

Consider  the  compressible  conservation  of  mass,  momentum  and  energy,  which  are  given 

by 


dp  d 
dt  dx{ 
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=  0 


dpuj 

dt 
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gXj  (  P<s«  +  T«) 

d 

[— Qi  +  (  P^ij  +  T~ij )  uj\ 


(la) 

(lb) 

(lc) 


22 


The  set  of  equations  given  by  Eqn.  (1)  is  closed  by  specifying  the  equation  of  state  p  =  pRT , 
where  R  is  the  gas  constant. 

For  the  hybrid  RANS-LES  formulation,  the  density,  pressure,  and  heat  transfer  coeffi¬ 
cients  in  the  Navier-Stokes  equations  remain  Reynolds- averaged,  while  the  remaining  vari¬ 
ables  are  Favre-  or  mass-averaged  to  include  compressibility  effects.  When  the  filtering 
operation  is  applied  to  the  Navier-Stokes  equations  (and  assuming  that  the  filter  commutes 
with  the  differentiation  operator,  which  is  strictly  not  correct),  the  general  form  of  the  LES 
equations  can  be  obtained.  Several  fluctuating  terms  arise  with  the  filtering,  indicating 
that  these  arc  the  turbulence  terms  that  must  be  blended  and  switched  between  URANS 
and  LES  in  different  regions  of  the  solution.  These  terms  include  the  viscous  stress  tensor, 
J)Tij  =  — pu”u ",  and  arising  from  the  decomposed  energy  equation,  the  turbulent  heat  flux, 


du'f 


3 1 


dxj' 


For  simulations  that  encom- 


qi  -  pu"h ",  and  the  rate  of  turbulent  dissipation,  pe  a 

pass  Mach  numbers  below  the  low  supersonic  regime,  the  molecular  diffusion  and  turbulent 
transport  terms  are  negligible  and  are  typically  not  included,  which  is  true  here. 

Closure  for  the  viscous  stress  tensor  is  obtained  from  the  Boussinesq  approximation  and 
is  modeled  using  the  typical  eddy  viscosity  approach 


2 

Txj  —  2  V’j'  S  ij  ~  (2) 

where  the  eddy  viscosity,  vt,  is  determined  from  the  URANS  and  LES  equations. 

The  Menter  SST  turbulence  model  solves  two  partial  differential  equations  describing  the 
turbulence  kinetic  energy  k  and  the  dissipation  rate  per  unit  turbulence  kinetic  energy  a>, 
which  is  then  related  to  the  turbulent  characteristic  length.  These  equations  are 
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(3a) 


(3b) 


where  the  RANS  superscript  indicates  a  variable  based  on  the  Reynolds-averaged  Navier- 
Stokes  formulation  that  will  be  blended  with  its  LES  SOS  counterpart  to  resolve  the  viscous 
stress  tensor  from  Eqn.  (2).  There  are  a  number  of  constants  in  Eqn.  (3),  which  are  defined 
as  /3*=0.09,  (7k  =0.85,  =0.5,  <7^2  =0.856,  and 


F2  =  tanh 
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where  y  is  the  distance  to  the  nearest  wall. 

The  turbulent  kinetic  energy  production  term 


k  is  computed  in  RANS  models  using 


either  vorticity  =  -y/fiyfi* 


where  =  - 
3  2 


1  (  dvii  du 


wThere  Sa  =  - 
13  2 


1  (  dui  duj 
dxj  dxi 


_ 2 

dxi 


,  or  the  strain  rate  S  =  yj2S~jSij 

y  l/  .  /  j  y 

Spalart70  noted  that  while  vorticity  has  been  more  commonly 


applied  in  many  RANS  models,  the  actual  physics  of  turbulence  production  are  based  on  the 
strain  rate.  In  the  boundary  layer  the  choice  of  vorticity  or  strain  rate  makes  little  difference, 
since  the  values  arc  relatively  comparable  in  that  region.  However,  in  other  regions,  such 
as  wakes,  where  vorticity  is  much  larger  than  the  strain  rate,  k  can  be  over  predicted  in  a 
vorticity-based  model.  Menter's  original  formulation69  in  1994  was  based  on  the  vorticity, 
but  his  more  recent  formulations71,72  recommend  a  change  to  the  strain  rate. 

Concurrent  with  the  solution  of  the  URANS  turbulence  model  equations,  a  one-equation 
LES  equation  for  the  turbulent  kinetic  energy  is  also  solved  to  determine  the  SGS  values  of 
the  blended  variables: 
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Here,  Ct  is  a  constant  0.916  and  A  is  the  characteristic  length  given  by  A  =  VJU  where  Vceu 
is  the  cell  volume  (for  the  unstructured  node-based  solver,  this  is  the  dual  cell  volume).  The 
SGS  eddy  viscosity  is  next  computed  from 


4GS  =  cv 


(6) 


where  Cv  is  a  constant  value  0.0667. 


5.1.4  Post  Processing  and  Analysis 

Post  processing  of  the  data  obtained  during  the  different  simulations  utilize  a  number  of 
different  software.  They  include  Tccplot73  and  Ficldvicw74  for  dynamic  and  static  views  of 
the  flow  field  and  surface  behavior,  and  MATLAB75  to  extract  and  plot  frequency  data. 
Data  and  contour  plots  are  presented  throughout  this  report.  In  these  plots,  the  following 
variable  definitions  were  used: 


Pressure  Coefficient:  Cp 

Q-eriterion:  Qcrit 

Vorticity  Magnitude:  \u\ 


V  ~  Poo 
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Analysis  of  the  computational  data  was  performed  to  ensure  that  the  minimum  error 
resulting  from  numerical  modeling  is  included  in  the  data.  In  section  5.2.4,  validation  of  the 
code  modifications  is  discussed,  along  with  verification  of  the  simulations.  Grid  independence 
studies,  within  the  scope  of  the  changing  scales  of  turbulence  eaptured/modeled  by  the  LES 
method,  have  been  undertaken.  Additional  levels  of  grid  adaptation  on  select  cases  were 
applied  to  discern  their  influence  on  performance  and  wake  characteristics. 

Statistical  data  analyses  were  also  performed.  Since  the  computational  results  cannot 
be  performed  for  the  large  number  of  revolutions  found  in  experiments,  a  moving  window 
analysis  was  made  on  averaged  data.  Windows  on  the  ultimate  and  penultimate  revolutions 
were  compared  to  ensure  that  the  solution  variable  of  interest  was  periodic  or  not  changing; 
if  not  additional  revolutions  were  simulated  until  this  behavior  was  observed.  The  results 
were  then  averaged  over  a  minimum  of  two  revolutions  to  provide  both  means  and  unsteady 
bounds1.  Care  was  taken  in  the  analysis  of  periodic  data  to  include  a  full  cycle  during  the 
analysis  to  eliminate  possible  bias  resulting  from  averaging. 

For  statistics  resulting  from  a  full  revolution,  an  additional  analysis  for  the  Gaussian  or 
normal  distribution  was  computed  to  determine  the  repeatability  of  the  unsteady  quantities. 
Data  lying  outside  the  tenth  and  ninetieth  percentiles  were  deemed  as  outliers.  Details  of 
these  analyses  are  provided  with  the  appropriate  data  results.  The  following  definitions  to 
compute  the  statistical  quantities  were  utilized,  with  the  MATLAB  version  R2013b  function 
name  included  m  parentheses: 


Mean  (mean): 

Standard  deviation  (std): 

Kurtosis  (kurtosis): 
Sample  Skewness  (skewness): 
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(11) 

(12) 

(13) 

(14) 


Fast  Fourier  Transforms  (FFT)  were  also  computed  to  find  the  turbulence  spectra  of  the 
data  at  specific  points.  For  consistency,  the  FFT’s  were  computed  using  the  MATLAB75 
function  fft.  This  function  (MATLAB  version  R2013b),  uses  the  formulae  based  on  the 
discrete  FFT  algorithm  by 


*  (*0  =  (15) 

;= i 

:An  unsteady  bound  should  not  be  confused  with  an  uncertainty  or  error  bound.  In  highly  unsteady 
flows,  such  as  those  encountered  here,  fluctuations  in  parameters  arc  expected,  and  their  quantification  is 
useful  in  many  engineering  applications. 
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where 


(16) 


fc=l 

and  where  ton  =  e“27r4/n. 


5.2  Development  of  Overset  Feature-Based  Grid  Adaptation 

To  help  meet  the  objectives  and  goals  of  the  computational  portion  of  this  research,  a  new 
adaptation  strategy  has  been  developed  that  permits  time-dependent  anisotropic  adaptation 
for  dynamic  overset  simulations.  This  strategy  minimizes  numerical  dissipation  and  compu¬ 
tational  costs  due  to  prohibitively  large  grid  sizes.  A  clear  advantage  of  this  technique  is 
that  grids  may  be  optimized  for  different  flow  conditions  and  geometric  orientations  with¬ 
out  need  to  manually  tailor  the  grid  each  time  a  new  run  condition  is  simulated.  Previous 
numerical  investigations5  76  of  hub  drag  have  not  examined  the  ability  of  grid  adaptation  to 
improve  predictions  of  the  loading  or  the  flow  features  of  the  unsteady  wake,  and  this  is  a 
new  capability  of  direct  use  to  the  Navy. 

The  current  development  permits  adaptation  to  be  executed  over  a  periodic  time  win¬ 
dow  in  a  dynamic  flow  field  so  that  an  accurate  evolution  of  the  unsteady  wake  may  be 
obtained  within  a  single  unstructured  methodology.  Unlike  prior  adaptive  schemes,  this 
approach  permits  grid  adaptation  to  occur  seamlessly  across  any  number  of  grids  that  arc 
overset,  excluding  only  the  boundary  layer  to  avoid  surface  manipulations.  Using  a  rotor- 
fuselage  interaction,  flow  field  physics,  time-averaged  and  instantaneous  fuselage  pressures, 
and  wake  trajectories  have  been  computed  and  compared  with  experiment.77  The  ability  of 
the  methodology  to  improve  these  predictions  without  user  intervention  has  been  confirmed 
(see  Section  5.2.4). 

The  ability  to  perform  adaptive  mesh  h—refincrrient  on  a  single  grid  is  made  possible  by 
the  refine  library.  Adaptive  mesh  refinement  (AMR)  can  be  performed  via  serial  or  parallel 
execution The  application  of  this  library,  prior  to  this  effort,  has  included  sonic-boom 
propagation,40  41  viscous  transonic  drag  prediction,40  and  re-entry  vehicle  configurations.42 


5.2.1  Anisotropic  Feature-Based  Adaptation 

FUNSD’s  anisotropic  tetrahedral  adaptation  capability39, 78  forms  the  basis  for  the  new  adap¬ 
tation  strategy.  This  feature-based  adaptation  requires  the  identification  of  a  feature,  as  well 
as  an  algorithm  or  key  to  define  the  grid  modification.  In  this  effort,  the  indicators  explored 
were  vorticity,  pressure  difference,  and  the  Q-critcrion.  The  vorticity  adaptation  formula¬ 
tion,  Fe^ |,  is  similar  to  the  one  applied  by  Duque  et  al., 38  which  scales  vorticity  (w)  with 
an  edge  length,  le.  For  a  given  edge  connecting  nodes  n\  and  n 2,  the  vorticity  formulation 
is  computed  based  on  the  averaged  vorticity  magnitude  across  the  edge, 


+  Mn2 

2 


(17) 
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The  pressure  difference  formula,  Fe^p,  was  defined  as  the  magnitude  of  the  pressure  difference 
A p  over  an  edge  scaled  by  the  edge  length  (le)  as 


^e,A p  le\Pnx  Pri2 1  ■  (1®) 

The  formulation  of  the  Q-criterion  indicator  is  based  on  Kamkar’s  non-dimensional  method79 
and  uses  the  rotation  rate  (f2)  and  strain  rate  (S)  tensors.  Here  the  maximum  value  across 
the  edge  is  applied, 

-SS^W"1))-  '  (19) 

Using  one  of  these  formulations,  the  normalized  local  adaptation  intensity,  7,  is  derived  for 
each  node  as  the  maximum  of  the  edge  key,  Ke,  over  all  incident  edges  of  a  given  node, 

/  =  max  ( -jr\ ,  (20) 

edges  y  / 

where  Kt  is  a  user-specified  tolerance.  The  new  isotropic  mesh  spacing  is  calculated  using 
an  estimate  of  the  spacing  on  the  original  mesh  h°,  a  coarsening  factor  C  (typically  around 
115%),  and  the  adaptation  intensity  by, 


hi  =  h°  min 


(21) 


The  power  of  0.2  is  an  under- relaxation  parameter  that  controls  the  aggressiveness  of  the 
refinement  process.  This  parameter  relates  the  convergence  rate  of  error  to  the  grid  spacing 
for  adjoint-based  adaptation.  The  value  of  0.2  (or  1/5)  has  been  found  sufficient  since  the 
convergence  rates  for  adjoint-based  adaptation  are  about  0(hA)  —  0(/i5).80  Although  there 
is  no  formal  connection  to  the  local  error  estimates,  this  value  has  been  typically  chosen  as 
a  sufficient  under-relaxation  parameter  for  feature-based  adaptation 42  utilized  in  this  effort. 

Consequently,  an  anisotropic  adaptation  metric  may  be  derived  using  a  scalar  quantity 
for  the  isotropic  spacing  and  a  Hessian  to  stretch  the  resulting  mesh.  The  Hessian  of  a 
quantity  (— )  can  be  described  as 


(22) 


Further  details  of  the  computation  Hessian-based  metric  and  its  significance  to  the  adapta¬ 
tion  process  are  delineated  in  Ref.  78.  The  vorticity-based  adaptation  invokes  the  vorticity- 
magnitude  Hessian  to  determine  anisotropy,  while  the  pressure  gradient  adaptation  utilizes 
the  Mach  number  Hessian,  described  in  Ref.  42. 
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5.2.2  Extension  to  Overset  Grids 

The  overset  grid  adaptation  capability  docs  not  restrict  adaptation  to  any  component  grid. 
This  enables  each  grid  to  evolve  independently  and,  in  general,  ensures  for  an  orphan  free 
composite  grid.  The  current  adaptation  capability  for  viscous  flows  is  restricted  only  to  nodes 
beyond  the  boundary  layer.  FUN3D  applies  an  adaptation  software  module  that  computes 
the  adaptation  metric  as  well  as  handles  all  associated  adaptive  meehanies.  Description  of 
the  parallelized  adaptation  mechanics,  which  include  grid  operations  such  as  node  insertion 
and  removal  by  splitting  or  collapsing  edges,  edge  and  face  swapping,  and  node  smoothing 
are  detailed  in  Ref.  78.  The  extension  of  this  method  to  include  overset  adaptation  requires 
communication  with  DiRTlib,59  the  grid  connectivity  module,  to  assign  a  component  mesh 
ID  for  each  node  in  the  composite  mesh.  The  code  performs  adaptation  over  the  entire 
composite  grid  system  by  tracking  the  component  mesh  ID  for  all  added  nodes. 

Since  overset  assembly  of  the  component  meshes  is  handled  by  a  library  outside  of  the 
FUN3D  framework  (SUGGAR++),49  a  generalized  global  index  convention  was  requisite 
so  that  subsequent  assembly  of  the  adapted  grid  with  its  domain  connectivity  information 
would  be  compatible  with  the  solution  information.  This  process  is  required  to  perform 
valid  solution  transfers  between  the  unadapted  and  adapted  grid  systems.  The  convention 
requires  both  the  flow  solver  and  adaptation  code  to  assign  composite  grid  global  indexes  by 
arranging  nodes  in  contiguous  fashion  by  mesh  ID  over  the  list  of  component  meshes.  Nodes 
added  due  to  adaptation  arc  initially  assigned  new  global  indexes  by  appending  them  to  the 
current  global  index  list.  Node  removal  results  in  unused  global  indexes,  which  is  handled 
by  a  reverse,  global-index  shifting  procedure.  In  order  to  satisfy  the  condition  of  contiguous 
mesh  IDs,  a  new  procedure  was  introduced  to  re-sort  the  global  indexes  of  the  adapted  grid 
system  as  illustrated  in  Fig.  2.  After  adaptation,  the  component  meshes  are  then  saved, 
and  the  resultant  domain  connectivity  information  is  obtained  by  invoking  SUGGAR++  for 
subsequent  grid  assembly. 
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(a)  Node  removal  (struck  through)  and  insertion  (b)  Resorted  global  indexes 
Figure  2:  Global  index  convention  illustrated  for  an  example  overset  grid  system. 


5.2.3  Time-Dependent  Adaptation 

Time-dependent  adaptation  is  obtained  using  a  methodology  based  on  that  developed  by 
Alauzet  and  Olivier.47  The  anistropic  grid  metric  is  computed  for  every  grid  node  at  a 
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given  time  step  and  is  progressively  intersected  over  a  selected  time  window  such  that  the 
strongest  restrictive  metric  at  each  node  is  retained  to  form  the  time-dependent  grid  metric. 
The  Hessian,  H ,  is  intersected  in  time  by  collecting  N  solution  samples  within  a  time  window 
w  as  given  by, 

N 

\Mw.max\  =  P|  |M*>n|.  (23) 

n= 1 

Using  the  resulting  metric,  a  new  adapted  mesh  may  be  obtained  suitable  for  multiple  time 
intervals  characterized  by  the  flow  phenomena  obtained  within  the  adaptive  window.  For 
rigid-body  rotating  systems,  such  a  window  can  be  identified  as  the  time  corresponding 
l/ribiades  revolutions.  The  solution  is  sampled  at  each  time  step  throughout  the  time  window 
to  obtain  the  metric  intersection.  The  present  approach  does  not  utilize  a  solution  transfer 
capability;  hence,  no  additional  interpolation  errors  are  introduced  and  the  inherent  solver 
accuracy  is  retained  as  a  result  of  re-simulation  of  the  unsteady  physics. 

5.2.4  Method  Validation 

The  new  overset  anisotropic  grid  adaptation  approach  was  validated  using  two  eases:  a 
rotor- fuselage  interaction  and  a  rotor  wake. 

Rotor-Fuselage  Interaction 

The  Georgia  Institute  of  Technology  (GIT)  rotor-fuselage  interaction  (RFI)  configuration 
is  comprised  by  a  cylindrical  fuselage  and  a  hemispherical  nose  to  permit  easier  identification 
of  RFI  has  been  extensively  evaluated  in  the  Harper  Wind  Tunnel.81  The  rotor  blades  have 
a  rectangular  planform  with  a  NACA-0015  airfoil  section.  The  rotor  blades  arc  nearly 
rigid  which  allow  for  CFD  analyses  that  neglect  structural  deformations.  Two  advance 
ratios  (//,  =  U.10  and  ft  =  0.20)  are  selected  for  investigation  and  the  relevant  blade  angles 
and  thrust  arc  reported  in  Tabic  2.  Data  from  this  effort  include  instantaneous  and  time- 
averaged  pressures  along  the  fuselage,  as  well  as  vortex  behavior  via  laser  light  sheets.  The 
fuselage  length  is  non-dimensionalized  (x/R)  by  the  rotor  radius  (7?  457mm)  for  ease  in 

presentation. 

This  model  has  been  evaluated  by  numerous  prior  computational  efforts  with  a  variety  of 
approaches,  for  example  Refs.  22,23,82-84.  O’Brien24  used  this  as  a  validation  ease  for  his 
series  of  actuator  to  overset  rotor  models  implemented  into  FUN3D.  Numerous  details  of  the 
time- aver  aged  fuselage  pressure  coefficient  have  not  been  captured  by  these  methods,  in  spite 
of  the  simplistic  model  geometry.  O’Brien24  noted  that  some  time-averaged  features  just  aft 
of  the  rotor  were  captured  when  the  entire  model  (rotor  strut  and  hub)  were  included,  which 
prior  efforts  neglected.  However,  the  vortex  interaction  observed  in  the  original  experiments 
near  the  nose  (x/R  0.3)81  has  not  been  adequately  resolved  by  any  of  these  prior  simulations. 

Simulations  were  computed  using  the  compressible,  inviscid  equation  set  as  separation 
and  other  viscous  effects  should  be  minimal  for  the  configurations  chosen.  Solution  ad¬ 
vancement  was  performed  with  a  time  step  equivalent  to  1°  azimuthal  sweep.  During  each 
time-step,  25  subiterations  were  used  in  conjunction  with  the  temporal  error  control  option 
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Pis 

Pis 

Ct 

0.10 

-2.02 

-1.94 

0.009045 

0.20 

-2.62 

-3.29 

0.009950 

Tabic  2:  Blade  angles  and  resulting  thrust  obtained  for  GIT  test  eases.  (Note:  All  angles 
are  reported  in  degrees.) 

to  ensure  one  to  two  orders  of  magnitude  reduction  in  residual.  The  metric  intersection 
was  performed  over  a  time  window  corresponding  to  180°  blade  sweep  or  180  steps  after  the 
solution  became  periodic  (after  two  revolutions).  In  order  to  obtain  valid  comparisons  with 
experiment,  the  thrust  values  obtained  for  these  cases  were  ascertained  to  be  within  1%  of 
those  listed  in  Table  2. 

A  grid  independence  study  on  the  fuselage  mesh  revealed  very  little  variation  in  steady 
surface  pressures.  This  evaluation  was  performed  using  the  same  near-body  mesh  for  the 
rotor  blades,  consisting  of  2.3  million  nodes  (13.5  million  cells).  The  baseline  fuselage  (back¬ 
ground)  mesh  consisted  of  5.1  million  nodes  (30.2  million  cells)  -  coarser  meshes  were  not 
considered  since  grid  resolution  around  the  rotor-disk  region  needed  to  be  maintained  for 
spacing  compatibility  with  the  blade  meshes.  Steady  state  fuselage  pressures  on  a  refined 
fuselage  mesh  consisting  of  9.6  million  nodes  (56.7  million  cells)  exhibited  minor  variation 
from  the  pressures  of  the  baseline  mesh.  Grid  convergence  for  the  baseline  mesh  was,  thus, 
established  and  its  application  for  the  full  rotor- fuselage  simulations  was  substantiated. 


Figure  3:  Model  of  the  GIT  rotor-airframe  configuration. 

Several  feature-based  adaptations  were  evaluated  for  the  GIT  RFI  model  to  determine 
the  validity  of  the  method,  as  well  as  the  appropriate  flow  field  metric.  The  most  accurate 
metric  was  determined  to  be  a  combination  of  the  pressure  gradient  and  vorticity  for  inviscid 
simulations.  Details  of  this  study  can  be  found  in  Shenoy  ct  al.77 

The  capability  of  this  adaptation  methodology  to  capture  viscous  flow  phenomena  was 
considered  since  this  is  important  for  the  hub  study  in  this  effort.  Two  turbulence  meth¬ 
ods,  the  URANS  model  (kcj-SST)  and  the  HRLES  model,  were  studied.  Using  a  highly 
pre- refined  composite  grid  (15.4  million  nodes),  where  the  background  grid  was  refined  in 
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the  wake  region  between  the  blade  and  the  fuselage,  the  effect  of  tip  vortex  dissipation  was 
studied  by  applying  both  these  turbulence  models.  The  eddy  viscosity  and  vorticity  mag¬ 
nitude  predicted  using  the  kw-SST  model  and  HRLES  models  arc  showm  in  Figs.  4  and  5, 
respectively.  These  figures  depict  the  prediction  at  ip  =  120°,  but  the  same  result  is  obtained 
at  different  blade  azimuths. 


l - x  l - X 


(a)  Eddy  viscosity  (b)  Vorticity  magnitude 

Figure  4:  Contours  along  the  airframe  symmetry  plane  from  the  kcj-SST  model. 
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(a)  Eddy  viscosity  (b)  Vorticity  magnitude 


Figure  5:  Contours  along  the  airframe  symmetry  plane  from  the  IIRLES  model. 

Two  observations  can  be  made  with  respect  to  Figs.  4  and  5.  First,  the  eddy  viscosity 
prediction  from  the  ku;-SST  simulation  is  significantly  higher  and  widespread  in  the  region 
coinciding  with  the  forward  tip  vortex.  The  vortex  also  appears  to  be  visibly  diffused  or 
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spread  out  in  comparison  to  the  HRLES  simulation.  Furthermore,  the  vortex  core  region  of 
the  HRLES  simulation  predicts  significantly  lower  eddy  viscosity,  implying  that  this  vortex 
core  exhibits  expected  laminar  behavior.  The  other  observation  is  that  the  high  eddy  vis¬ 
cosity  prediction  from  the  ko;-SST  simulation  that  dominates  the  rotor  wake  region  clearly 
diffuses  the  vortex  propagating  toward  the  fuselage  following  interaction  with  the  oncoming 
blade.  Additionally,  rotor  wake  vortieity  contours  of  the  HRLES  simulation  show  two  dis¬ 
tinct  high  vortieity  regions  that  are  characteristic  of  the  expected  vortex-fuselage  interaction. 
This  study  demonstrates  the  superiority  of  the  HRLES  model  in  the  preservation  of  the  tip 
vortex  and  capturing  a  more  complex  interaction  with  the  subsequent  blade  passage. 

The  initial  grid  used  for  the  HRLES  computations  had  the  same  pre-refinement  in  the 
rotor-wake  region  as  in  the  inviseid  simulations.  The  grid  sizes  resulting  from  the  initial  grid 
to  the  different  feature-based  schemes  are  listed  in  Table  3.  Among  the  single  adaptation 
schemes,  adaptation  to  A p  was  not  considered  because  of  its  inability  to  preserve  regions  of 
high  vortieity,  which  are  essential  in  order  to  capture  the  magnitude  of  the  fuselage  surface 
pressures.  The  double  adaptation  schemes  studied  were  the  vorticity-mixed  scheme  and  the 
Q- criterion-mixed  scheme.  The  |o;|  (iter.  2)  scheme  was  not  performed  because  the  inviseid 
simulations  showed  very  little  improvement  from  the  single  adaptation  to  \u\. 

Table  3:  HRLES:  summary  of  grid  size  resulting  from  the  different  adaptation  schemes 


Scheme  Description 

Total  Nodes  (Millions) 

Initial  Grid 

5.52 

Adapted  to  |cj| 

13.0 

Adapted  to  Q-erit. 

16.8 

Adapted  to  |cj|  &  A p 

24.2 

Adapted  to  Q-erit.  &  A p 

30.7 

Time-averaged  fuselage  pressures  comparing  the  \uj\  and  Q-eriterion  single  adaptation 
schemes  (Fig.  6)  show  very  small  differences,  unlike  the  results  from  the  inviseid  simulation.77 
Both  these  schemes  improve  the  magnitude  of  the  pressures  in  the  forward  part  of  the 
fuselage  (x/R  <  0.3),  with  exception  to  the  nose  region.  Additionally,  the  primary  vortex 
interaction  resulting  from  both  these  schemes  indicate  the  presence  of  a  small  pressure  pulse 
at  x/R  «  0.5. 

The  vorticity-mixed  scheme  and  Q- criterion-mixed  scheme  improve  the  accuracy  of  the 
time-averaged  pressures  (Fig.  7),  particularly  in  predicting  the  intensity  of  the  primary  vortex 
interaction  x/R  «  0.5.  Both  schemes  have  clear  similarities,  but  the  vorticity-mixed  scheme 
predicts  a  stronger  pressure  pulse  than  the  Q- criterion- mixed  scheme.  Again,  since  the  hub 
pin  geometry  is  not  modeled,85  the  HRLES  simulations  are  not  able  to  accurately  predict 
the  pressures  in  the  aft  portion  of  the  rotor  (x/R  >  1.5). 

The  instantaneous  fuselage  pressures,  plotted  in  Figs.  8  and  9,  show  good  correlation  with 
experimental  data  at  intermediate  azimuths  (ij)  =  90°  —  150°)  with  both  the  vorticity-mixed 
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x  R 


Figure  6:  IIRLES:  time-averaged  pressures  on  the  top  centerline  from  one  adaptive  cycle. 


Figure  7:  IIRLES:  time-averaged  pressures  on  the  top  centerline  from  two  adaptive  cycles. 


scheme  and  Q- criterion-mixed,  scheme.  Both  schemes  agree  very  well  with  each  other  and 
the  only  difference  observed  is  that  the  vorticity -mixed  scheme  predicts  a  higher  magnitude 
of  the  primary  vortex  interaction.  With  respect  to  experiment,  both  schemes  result  in  a 
small  lead  (A x/R  =  0.01)  in  the  spatial  location  for  the  primary  vortex  interaction.  The 
magnitude  of  the  interaction  is  under-predicted  for  the  first  quarter  revolution  and  is  over- 
predicted  for  the  second  quarter  revolution.  The  secondary  vortex  interaction  initially  shows 
a  slight  lead  in  the  spatial  location  (Ax/ R  =  0.05),  but  during  the  second  quarter  revolution, 
this  interaction  lags  behind  the  experiment  by  the  same  amount.  The  faster  convection  rate 
of  the  secondary  vortex  interaction  is  not  captured.  This  may  be  attributable  to  boundary 
layer  effects  not  being  modeled  as  accurately  either  due  to  lack  of  adaptation  or  shortcomings 
of  the  turbulence  modeling  in  this  region.  As  with  the  inviscid  simulations,  the  high  pressure 
regions  in  the  forward  portion  of  the  fuselage  (0.2  <  x/R  <  0.5)  at  ip  =  180°  do  not  correlate 
well. 
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Figure  8:  HRLES:  top  centerline  instantaneous  pressures  (first-quarter  revolution). 
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Figure  9:  HRLES:  top  centerline  instantaneous  pressures  (second-quarter  revolution). 


To  further  understand  the  significance  of  the  grid  refinement  in  each  scheme,  the  vortex 
behavior  is  examined  in  Figs.  10  -  12.  The  adaptation  sequence  from  the  initial  grid  to  the 
vorticity-mixed  scheme  arc  plotted  from  top  to  bottom  at  each  selected  azimuthal  location. 
Since  the  adaptation  sequence  leading  to  Q-cntenon-mixed  scheme  show  the  similar  results 
as  with  the  vorticity-mixed  scheme ,  they  are  not  presented  here.  It  is  clear  from  scanning 
from  top  to  bottom  that  the  forward  vortex  core  is  more  crisply  predicted  after  the  first 
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adaptation  (middle  plot)  and  refined  further  upon  the  second  adaptation  (bottom  plot).  In 
addition,  the  initially  weaker  vortical  features  in  the  rotor  wake,  diffused  over  large  areas, 
become  further  defined  due  to  adaptation. 


x/R  x;R 

(a)  ip  =  210°  (ip  =  30°)  (b)  ip  =  240°  (ip  =  60°) 

Figure  10:  Vortex  behavior  (from  top:  initial  grid,  adapted  to  |o;|,  and  vorticity-mixed 
scheme ). 


36 


■ 


o.2  o.4  oo  o.e 

x/n 
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(b)  ip  =  300°  (ip  =  120°) 


Figure  11:  Vortex  behavior  (from  top:  initial  grid,  adapted  to  |(u|,  and  vorticity-mixed 
scheme ). 
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(a)  ip  -  330°  (ip  =  150°)  (b)  ip  =  360°  (ip  =  180°) 

Figure  12:  Vortex  behavior  (from  top:  initial  grid,  adapted  to  |u/|,  and  vorticity-mixed 
scheme). 

Differences  between  the  various  adaptation  schemes  can  also  be  discerned  from  the  mag¬ 
nitude  and  shape  of  the  vortieity  contours  in  these  figures.  For  example,  the  shape  of  the 
vortex  at  0.4  <  x/R  <  0.5  where  the  primary  interaction  occurs  is  very  different.  Specifically, 
tracing  the  vorticity-mixed  scheme  across  Figs.  10  -  12,  it  is  possible  to  discern  the  path 
of  the  tip  vortex  as  it  leaves  the  blade,  interacts  with  the  previous  blade’s  wake  sheet,  and 
finally  collides  with  and  encompasses  the  fuselage  centerline.  Brand81  reported  that  the  tip 
vortex  from  the  prior  blade  interacts  with  the  following  blade  at  x/R  =  0.3  at  ip  =  188°  (or 
ip  =  8°),  which  is  comparable  to  the  grid  adaptation  results  in  Fig.  12  (b).  The  weakening 
vortieity  of  the  primary  interaction,  observed  at  approximately  x/R  0.45  in  Figs.  10  -12, 
correlates  to  the  experimental  visualization.  The  vortex  sheet  roll-up,  which  was  experimen- 
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tally  observed  to  traverse  in  the  range  0.4  <  x/R  <  0.5  during  this  azimuthal  time  period 
can  also  been  observed  traveling  downstream  at  ijj  =  210°  and  with  a  distinct  rotation  by 
=  240°,  located  in  the  same  fuselage  locations.  The  development  of  the  secondary  vortex- 
fuselage  interaction  and  its  subsequent  rapid  downstream  convection  is  also  observed  in  the 
vorticity-mixed  scheme  and  Q- criterion-mixed  scheme  (not  shown). 

A  lack  of  agreement  between  the  vorticity-mixed  and  Q-critenon-mixed  schemes  was 
observed  in  inviseid  simulations77  but  not  in  the  viscous  simulations.  The  Q-eriterion  in¬ 
dicator  formulation  targets  regions  where  the  rotation  rate  \\Q\\  dominates  the  strain  rate 
||5j|,  since  Ftoi  =  0.01.  In  regions  where  H^H  exceeds  ||fi||,  the  Q-eriterion  (dimensional  or 
non-dimensional)  values  are  negative.  Therefore  such  regions  are  not  selected  for  refinement. 
The  vorticity  magnitude  method,  on  the  other  hand,  does  not  discriminate  regions  where 
strain  rates  dominate  and  its  range  is  always  non-negative. 

Figure  13  illustrates  the  differences  in  the  flow  fields  resulting  from  these  schemes  for 
the  inviseid  simulations.  The  vortex  sheet  region  clearly  displays  uniformly  high  values 
of  vorticity,  but  the  values  of  Q-eriterion  in  those  regions  are  negative  or  very  close  to 
zero.  Regions  where  both  vorticity-magnitude  and  Q-crit.erion  are  high  include  the  vortex 
core  region  and  few  localized  regions  in  the  vortex  sheet.  The  fully  turbulent  flow  fields 
resulting  from  these  adaptation  schemes  is  compared  in  Fig.  14.  High  values  of  vorticity 
throughout  the  vortex  sheet  are  observed  in  Fig.  14  (a).  In  comparison  to  Fig.  13  (b), 
the  vortex  sheet  in  Fig.  14  (b),  shows  several  more  spots  of  positive  Q-eriterion  values. 
Therefore,  these  regions  in  the  vortex  sheet  are  selected  for  refinement  in  the  turbulent 
simulation  and  are  excluded  from  refinement  in  the  inviseid  simulation.  Since  the  adaptation 
metric  is  accumulated  over  180°  azimuthal  sweep,  the  highly  refined  vortex  sheet  benefits 
the  preservation  of  the  tip  vortex  as  it  passes  through  these  regions.  Another  difference 
between  inviseid  and  turbulent  simulations  is  that  the  latter  exhibits  an  interaction  with 
the  rotor  wake  and  fuselage  boundary  layer  about  the  juncture  of  the  hemispherical  nose 
and  cylindrical  fuselage.  This  region  indicates  both  high  values  of  vorticity  magnitude  and 
positive  Q-eriterion.  This  behavior  arises  mainly  due  to  the  no-slip  boundary  condition  of 
the  turbulent  simulations. 
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(a)  Vo rticity- mixed  scheme  (vortirity-  (b)  Q- criterion-mixed  scheme  (Q- 
magnitudc  contours)  criterion  contours) 


Figure  13:  Inviscid  simulations:  comparison  of  the  vorticity-mixed  and  Q- criterion- mixed 
schemes. 


(a)  vorticity-mixed  scheme  (vortieity-  (b)  Q-cntenon-mixed  scheme  (Q- 
magnitude  contours)  criterion  contours) 


Figure  14:  HRLES:  comparison  of  the  vorticity-mixed  and  Q- criterion- mixed  schemes. 

The  genesis  of  positive  Q-criterion  values  in  the  turbulent  vortex  sheet  is  explained  by 
examining  the  flow  field  near  the  blade  trailing  edge  in  Fig.  15.  The  trailing  edge  region 
exhibits  high  Q-criterion  values  over  a  significant  portion  of  the  blade  span  for  the  turbulent 
simulation,  absent  from  the  inviscid  simulation.  This  is  attributed  to  boundary  layer-trailing 
edge  vortex  shedding  due  to  the  blade’s  no-slip  boundary  condition. 

Experimental  visualizations  document  the  tip  vortex  locations  as  soon  as  they  become 
visible  at  —  188°  and  are  plotted  at  30°  intervals  until  impingement  with  the  fuselage 
in  Fig.  16  (a).  Vortex  locations  from  the  vorticity-mixed  scheme  resulting  from  both  the 
inviscid  and  HRLES  simulations,  are  shown  in  comparison  to  experiment.  Figure  16  (b)  plots 
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(a)  Inviscid  Q- criterion-mixed  scheme  (b)  IIRLES  Q- criterion-mixed  scheme 


Figure  15:  Comparison  of  inviscid  and  viscous  modeling  of  the  flowfield  for  the  Q- criterion- 
mixed  scheme 


the  streamwise  spatial  location  lead  (or  lag)  with  respect  to  experiment.  Both  simulations 
show  the  same  lead  at  the  first  vortex  location.  However,  the  inviscid  simulation  shows  a 
vortex  lag  for  the  rest  of  the  azimuthal  locations,  as  high  as  Ax/ R  =  0.032.  The  IIRLES 
simulation  correlates  much  better  with  the  spatial  location,  generally  leading  the  experiment, 
with  the  maximum  vortex  lead  of  A x/R  —  0.013.  This  spatial  lead  is  also  observed  via  the 
surface  pressures  in  Figs.  8  and  9.  The  uncertainty  in  the  streamwise  location  reported 
with  experiment  was  15mm  or  A  x/R  =  0.033.  Both  simulations  with  the  vorticity-mixed 
scheme  arc  within  the  experimental  error,  but  Fig.  16  (b)  shows  significantly  better  spatial 
correlation  from  the  HRLES  simulation  of  the  vortex  location  with  experiment. 
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Figure  16:  Vortex-trajeetory  comparisons  from  the  vorticity-mixed  scheme. 


The  convergence  of  this  feature-based  adaptation  process  has  been  assessed  by  evaluating 
an  integrated  quantity  of  interest  obtained  over  a  series  of  adaptation  cycles.  Because  the 
time-averaged  pressure  distribution  has  been  previously  applied  as  a  suitability  criterion  of 
an  adaptation  scheme,  the  time-averaged  centerline  pressure  integral  is  used  here  to  identify 
convergence.  The  functional  quantity  is: 


(24) 


where  the  fuselage  length  is  3R. 

The  convergence  of  the  vorticity-mixed  scheme  has  been  assessed  by  performing  an  ad¬ 
ditional  adaptation  sequence,  i.e.  additional  adaptations  that  include  both  |u;|  and  A p.  The 
pressure  integral  is  plotted  in  Fig.  17  for  the  vorticity-mixed  HRLES  simulation.  Conver¬ 
gence  is  established  for  the  four  adaptation  cycles  since  the  change  in  the  functional  between 
cycles  3  and  4  is  observed  to  be  within  0.05%. 
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Figure  17:  Pressure-integral  functional  convergence  for  the  IfRLES  vorticity-mixed  scheme 
adaptation  sequence. 

The  efficacy  of  the  overset,  time- accurate  grid  adaptation  capability  has  been  demon¬ 
strated,  but  the  solution  dependency  on  the  selection  of  the  feature  or  combinations  of 
features  is  clear.  In  addition,  any  uncertainty  in  the  number  of  adaptation  cycles  suggests 
that  a  method  to  directly  relate  convergence  to  the  functional  of  interest  is  required,  such 
as  adioint-based  adaptation.  Adjoint-based  adaptation  has  been  demonstrated  for  steady, 
single  grids  by  Park  and  associated  authors.78, 86  The  extension  of  adjoint-based  adapta¬ 
tion  capability  to  include  overset,  time-accurate  simulations  may  be  warranted  and  is  being 
explored  beyond  this  project. 

Wake  Prediction 

As  this  effort  includes  the  analysis  of  rotor  hub  wakes,  a  second  validation  case  that 
includes  wake  data  was  evaluated.  The  second  validation  case  uses  the  ROtor  Body  IN- 
teraction  (ROBIN)  configuration,  developed  by  NASA,  which  has  been  extensively  used 
in  various  experiments  and  computational  studies  on  rotor-fuselage  interactions  and  wake 
trajectories.26  87-90  This  streamlined  slender  ROBIN  fuselage  model,  described  by  a  set  al¬ 
gebraic  equations  at  various  fuselage  locations,  yields  a  simple  analytical  definition  for  a 
fuselage  geometry.  An  engine  mount  or  “doghouse’’  is  a  characteristic  feature  included  in 
the  configuration.  An  internally  mounted  rotor  system  is  utilized  consisting  of  four  blades 
that  are  fully-articulated  with  a  NACA-0012  airfoil  section.  The  current  effort  focuses  on  the 
set  of  experiments  conducted  by  Ghee  and  Elliott88  in  the  14-by-22-Foot  Subsonic  Tunnel  at 
NASA  Langley  Research  Center  using  the  two-meter  rotor  test  system  (2MRTS).  One  of  the 
experimental  wake  visualization  cases  (Table  4)  was  employed  here  to  evaluate  the  influence 
of  the  adaptation  process.  Once  again,  the  length  data  have  been  nondimensionalized  by 
the  rotor  radius  (R)  to  facilitate  case  in  interpreting  the  simulations. 

The  compressible  viscous  option  within  FUN3D  was  applied  to  the  ROBIN  demonstra- 
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0.23 

0.0064 

-3.0 

1.5 

6.5 

-3.2  -1.1 

Table  4:  Relevant  parameters  for  the  ROBIN  wake  visualization  test  case.  (Note:  All  angles 
axe  reported  in  degrees.) 


tion.  Similar  to  the  GIT  computations,  the  time  step  equivalency  of  1°  was  chosen  with  25 
subiterations.  Ten  turbulence  subiterations  were  also  applied  to  converge  the  loosely-coupled 
Menter’s  kui-SST  turbulence  model.  Here  the  time  window  for  the  time-dependent  metric 
corresponded  to  90°  blade  sweep  or  90  steps.  Trim  angles  that  were  previously  obtained  on 
the  same  grid  in  Smith  et  al.26  were  used  to  estimate  the  required  trim.  The  thrust  values 
obtained  were  ascertained  to  be  within  2%  of  the  experimental  values  listed  in  Table  4,  which 
validate  the  correlations  with  experiment  without  further  trimming  the  rotor. 


Figure  18:  Computational  model  of  the  ROBIN  2MRTS  configuration. 

The  purpose  of  this  demonstration  case  was  to  examine  the  influence  of  grid  adaptation 
on  the  prediction  of  the  wake  vortex  trajectories.  To  examine  this,  a  baseline  grid  wake  tra¬ 
jectory  was  compared  to  a  single  adaptation  based  on  vorticitv  magnitude.  The  application 
of  only  a  single  adaptation  with  the  selection  of  the  vorticity  magnitude  as  the  adaptation 
metric  was  based  on  the  prior  GIT  RFI  simulations  for  the  range  of  experimental  data  avail¬ 
able  for  correlation.  The  initial  grid  had  14.4  M  nodes  and  the  vorticity  based  adaptation 
increased  the  size  to  17.6  M  nodes. 

Vorticity  tracks  were  extracted  at  different  longitudinal  slices  spanning  the  rotor  disk. 
Vortex  cores  are  determined  by  mapping  local  maxima  of  the  normal  vorticity  component  on 
each  of  these  planes.  There  was  no  attempt  made  to  differentiate  individual  vortex  structures 
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like  tip  vortices  from  these  contour  maps.  The  wake  trajectories  for  the  advance  ratio  of  0.23 
(Fig.  19)  are  observed  to  be  very  similar  to  the  experiment  and  the  baseline  grid  trajectories. 
The  major  influence  of  the  grid  adaptation  was  not  readily  transparent,  as  it  was  apparent 
primarily  in  the  process  through  which  the  wake  trajectories  were  numerically  determined. 
As  discussed  with  the  GIT  demonstration  case,  and  illustrated  here  via  Q-criterion  plots 
(Fig.  20),  the  effect  of  the  grid  adaptation  was  to  further  refine  the  vorticity  throughout  the 
flow  field.  Thus,  the  refined  strength  and  minimized  extent  of  the  vortex  path  rendered  the 
wake  tracking  much  simpler  and  with  fewer  approximations. 


0.5 


*  •  •  •  *  • 
* '  s'?  **  '  ■ 


•  Experiment 
■  Initial  Grid 

♦  Adapted  to  vort,  mag 


0,5 


0  0.5 

x/R 

(b)  y/R  =  0.3 


1.5 


(c)  y/R  —  —0.8 


Figure  19:  ROBIN  sideline  vortex  trajectories  at  //  —  0.23. 
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(b)  Adapted  Grid 


(a)  Initial  Grid 


Figure  20:  Q-critcrion  iso-surfaccs  colored  by  vorticity  magnitude  for  the  ROBIN  configura¬ 
tion  at  fi  =  0.23. 
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6  Hub  Analysis  and  Deconstruction 

6.1  Experiment  Description 

Experiments  on  a  generic  four-bladed  hub  model  were  conducted  in  the  John  J.  Harper  low 
speed  wind  tunnel  located  within  the  Daniel  Guggenheim  School  of  Aerospace  Engineering 
of  Georgia  Institute  of  Technology.  A  more  complete  detailed  description  with  analysis  can 
be  found  in  Part  1  of  this  report.52  The  model  (Fig.  21)  is  approximately  one-quarter  scale 
of  a  ten-ton  helicopter.  The  complete  model  includes  a  number  of  geometric  components 
found  in  a  typical  rotor  hub:  hub  plates,  blade  shanks,  swashplate,  pitch  links,  drive  shaft 
and  requisite  hardware  (nuts,  bolts,  etc).  A  modified  model,  tested  at  the  end  of  project 
in  2013,  also  added  a  scissors  assembly,  based  on  recommendations  from  industry.  Some 
components  found  in  full-scale  rotors,  such  as  hydraulic  lines  and  control  wire  bundles,  were 
not  included  in  the  model. 

While  there  were  different  models  with  plugged  and  unplugged  rotor  shanks,  only  the 
plugged  shank  model  was  tested  with  rotation.  A  third  configuration  included  plugged  rotor 
shanks  and  closed  the  open  region  between  the  two  hub  plates  (Fig.  21).  A  view  of  the 
experimental  model  and  mount  in  the  tunnel  test  section  is  illustrated  as  Fig.  22. 

The  experimental  evaluation  was  initiated  with  static  hub  tests.  Various  azimuthal  ori¬ 
entations  at  zero  angle  of  attack  were  tested.  Force  data  were  obtained  in  15°  increments 
over  a  quarter  revolution  (since  this  is  a  four-bladed  hub)  for  a  range  of  available  tunnel 
speeds.  The  plugged  shank  configuration  was  evaluated  over  the  similar  wind  tunnel  speed 
ranges  and  for  rotation  rates  up  to  240  rpm,  and  was  the  only  configuration  modeled  compu¬ 
tationally.  In  addition  to  the  force  measurements  on  the  model,  particle  image  veloeimetry 
(PIV)  and  singe  axis  hot-wire  probes  were  applied  at  selected  wake  locations  to  measure  the 
wake  velocity  and  frequency  spectra  of  velocity  fluctuations,  respectively. 

Hub  drag  deconstruction  was  achieved  experimentally  by  removing  various  components 
for  selected  static  and  rotating  conditions.  The  drag  was  measured  by  removing  the  blade 
shanks,  then  the  shanks  and  hub  plates,  pitch  links,  etc.  until  only  the  hub  drive  shaft 
remained. 
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(c)  Side  view  of  the  experimental  model  (Ref.  91) 


Figure  21:  Comparison  of  the  computational  (CFD)  and  experimental  plugged  geometries. 
The  CFI)  geometries  also  include  identification  of  components  and  orientation. 


48 


Figure  22:  Hub  model  within  John  J.  Harper  wand  tunnel  test  section.  From  Ref.  91. 


6.2  Grid  and  Run  Option  Descriptions 

Time-accurate  predictions  of  the  static  hub  and  unsteady  wake  were  simulated  using  the 
hybrid  Reynolds-Averaged  Navier-Stokes  and  large  eddy  simulation  (RANS/LES)  turbulence 
approach  (GT-HRLES)92  described  in  Section  5.1 .3.  The  background  grid  included  the  wind 
tunnel  test  section  (see  Fig.  22),  where  the  tunnel  walls  were  modeled  as  inviscid  surfaces. 
The  wind  tunnel  mount,  as  well  as  details  of  the  hub,  including  pitch  links,  bolts  and  nuts 
were  included  in  the  configuration,  as  previously  illustrated  in  Figs.  21  and  22.  Each 
component  was  identified  in  the  simulations  separately,  which  resulted  in  sets  of  integrated 
loads  for  each  component,  as  wTell  as  the  entire  assembly.  The  ability  to  group  and  monitor 
the  unsteady  history  of  loads  on  various  surfaces  is  beneficial  in  order  to  assess  the  drag 
predictions  against  theory  and  experiment. 

A  baseline  overset  grid  of  11.1  million  nodes  vras  adapted  using  vorticity  magnitude  as 
the  adaptation  indicator.  The  boundary7  layer,  consisting  of  35  normal  cell  layers,  had  a 
maximum  y+  —  0.35  for  the  reference  Reynolds  number  of  0.291  x  106.  An  example  of  a  grid 
before  and  after  adaptation  is  illustrated  in  Fig.  23.  The  importance  of  the  adaptation  is 
emphasized  in  the  shed  vorticity  of  Figs.  24  and  25.  The  adapted  grids  provide  significantly 
higher  fidelity  of  the  unsteady  shed  wake.  This  is  particularly  important  in  regions  where 
the  shed  wake  impacts  on  other  hub  components,  as  observed  for  the  shed  wake  of  the  hub 
upstream  of  the  main  strut  in  Fig.  25.  The  initial  and  final  adapted  grid  sizes  for  several  of 
the  cases  examined  computationally  are  listed  in  Table  5.  Rotating  hub  adaptations  resulted 
in  larger  overall  grids  than  their  static  counterparts  due  to  the  larger  extent  of  vorticity  in 
the  w^ake.  Typically  two  adaptation  cycles  were  required  to  result  in  drag  changes  of  less 
than  1%,  as  documented  by  Shenoy  et  al.93 


49 


In  addition  to  the  basic  computational  options  described  in  Section  5.1.2,  there  were 
several  options  that  pertains  specifically  to  the  hub  simulations.  A  dimensional  time  step 
(0.07  milliseconds)  was  selected  to  be  equivalent  to  a  one  degree  azimuthal  sweep  of  the 
240  rpm  rotating  hub.  Using  this  time  step,  a  fluid  particle  traverses  the  characteristic  hub 
length  in  about  80  time  steps.  Each  time  step  was  augmented  with  an  average  of  25  Newton 
subitcrations  (up  to  a  maximum  of  40  subitcrations)  to  increase  the  temporal  accuracy  of 
the  simulation.  A  temporal  error  controller  maintained  a  specific  residual  error  so  that  the 
number  of  subiterations  at  each  time  step  varies. 

To  mimic  the  turbulence  intensity  of  the  wind  tunnel,  a  free  stream  turbulence  intensity 
of  1%  was  applied  to  the  simulation.  This  value  was  chosen  based  on  prior-reported  tur¬ 
bulence  in  wind  tunnels  of  the  same  generation  as  the  Harper  wind  tunnel.  During  2012, 
there  were  major  renovations  to  the  wind  tunnel,  and  it  was  reported  that  the  wind  tun¬ 
nel  turbulence  intensity  reduced  to  0.03%.  However,  to  maintain  the  correlations  with  the 
pre-renovation  experiments,  the  original  turbulence  intensity  was  maintained.  Simulations 
on  similar  complex  configurations  with  both  turbulence  intensity  levels  did  not  result  in 
observable  changes  in  the  parameters  evaluated  in  this  report. 

The  free  stream  density  was  assumed  to  be  1.1901  kg/m3  with  a  reference  velocity  of 
8.941  m/s.  With  a  shank-to-shank  hub  diameter  of  0.4862  m,  the  Reynolds  number  was 
computed  to  be  0.291  x  106. 


(a)  Baseline  (Unadaptcd)  (b)  Anisotropic  grid  adaption 

Figure  23:  Application  of  anisotropic  grid  adaptation  for  increased  resolution  of  wake  field 
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(a)  Initial  rotating  grid  (b)  Adapted  rotating  grid 

Figure  24:  Planform  (looking  down)  view  of  the  vorticity  magnitude  for  the  initial  and 
adapted  grid. 


(a)  Initial  rotating  grid 


(b)  Adapted  rotating  grid 


Figure  25:  Left  (looking  forward)  side  view  of  the  vorticity  magnitude  for  the  initial  and 
adapted  grid. 
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Test  condition 

Initial  Grid 

Adapted 

Static  hub  at  0°  orientation  (U^  =  8.941  m/s) 
Static  hub  at  45°  orientation  (U^  =  8.941  m/s) 
Rotating  hub  at  240  rpm  (U^  =  8.941  m/s) 
Static  hub  at  0°  orientation  ( U, ^  =  13.41  m/s) 
Rotating  hub  at  240  rpm  (f/oo  =  13.41  m/s) 
Rotating  hub  at  240  rpm  (U^  =  22.35  m/s) 

11.1  M  nodes 
11.1  M  nodes 
11.1  M  nodes 
11.1  M  nodes 
11.1  M  nodes 
11.1  M  nodes 

15.8  M  nodes 

16.0  M  nodes 

19.5  M  nodes 

20.9  M  nodes  (2  iterations) 
29.0  M  nodes  (2  iterations) 

28.9  M  nodes  (2  iterations) 

Table  5:  Effect  of  grid  adaptation  on  the  grid  size. 

6.3  Error  Analysis 

With  any  investigation  there  are  errors  associated  with  both  experimental  and  computa¬ 
tional  efforts.  Errors  should  be  quantified  and,  if  possible,  minimized  prior  to  analysis  of 
data.  An  assessment  of  both  the  experimental  and  computational  errors  are  addressed  to 
provide  a  reference  point  for  these  and  future  analyses. 

6.3.1  Experimental 

Experimental  errors  associated  with  the  apparatus  have  been  defined  in  the  experimental 
portion  of  this  report.52  In  addition  to  the  documented  vendor  errors,  there  are  other  errors 
associated  with  the  experiment.  These  are  primarily  based  on  human-related  errors  that 
can  influence  the  measurements.  These  errors  may  include,  but  not  limited  to,  issues  with 
calibration  of  the  equipment  before  and  after  the  test,  measurement  of  the  locations  where 
the  wake  measurements  were  taken,  data  reduction  errors  (associated  with  approximations), 
and  setting/transeription  of  the  experimental  operational  conditions  (temperature,  pressure, 
velocity,  and  their  fluctuations  during  the  test).  Their  estimates  of  some  of  the  errors  are 
shown  in  Table  6. 


Parameter 

Error  Estimate 

Dynamic  Pressure 

0.02% 

Dvnamie  Viscosity 

0.01% 

Density 

0.03% 

Velocity 

0.03% 

Load  Cell  Calibration 

1.37% 

Calibration 

0.77% 

Temperature  variation 

0,45% 

Approximation 

0.43% 

Ambient  pressure 

0.03% 

Table  6:  Estimation  of  the  error  for  some  experimental  parameters.  From  Ref.  52 
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6.3.2  Computational 

Numerical  errors  fall  into  two  major  categories:  truncation  errors  and  floating  point  precision 
errors.  The  former  are  primarily  defined  in  the  development  of  the  computational  technique, 
although  they  can  be  minimized  with  a  refined  isotropic  grid.  The  second  is  associated  with 
the  computation  itself. 

The  use  of  the  feature-based  overset  grid  adaptation  aids  to  minimize  these  errors  as 
the  spatial  grid  is  refined  in  the  area  where  the  maximum  changes  (and  hence  the  largest 
gradients)  occur  in  the  flow  field.  To  minimize  these  errors,  grid  adaptation  was  performed 
until  the  primary  parameters  of  interest  (hub  drag,  as  well  as  velocities  and  turbulent  spectra 
at  specified  wake  locations)  do  not  change  within  a  defined  tolerance.  Examples  with  the 
hub  drag  error  minimization  arc  provided  as  part  of  this  report.  Since  the  turbulent  wake  is 
modeled  wfith  the  large  eddy  simulation  approach,  grid  refinement  will  result  in  the  capture  of 
more  of  the  turbulence  scales.  Thus,  changes  in  the  wake  behavior  will  change  until  the  grid 
has  been  minimized  such  that  all  of  the  turbulent  scales  that  contribute  to  the  parameter 
evaluated  have  been  captured.  This  is  a  direct  relation  to  the  tolerance  chosen  to  define 
convergence. 

Associated  with  the  temporal  integration  arc  both  truncation  errors  and  floating  point 
precision  errors.  In  order  to  minimize  computational  times,  most  unsteady  computational 
methods,  including  FUN3D,  apply  implicit  time  integration  (specifically  here  the  backward 
difference  formulation).  Convergence  of  the  time  step  is  demonstrated  in  Fig.  26,  where  the 
time  steps  studied  corresponded  to  0.5°  (At  ~  0.3472  x  10~3  s)  and  1°  (At  =  0.6944  x  10-3 
s)  azimuthal  sweeps.  Prediction  of  the  4/rev  (16  Hz)  and  8/rev  (32  Hz)  features  that  are 
prominent  at  these  wake  locations  is  insensitive  to  time  step  refinement  (Ref.  91).  Further, 
for  frequencies  higher  than  100  Hz,  a  broad  spectrum  of  scales  is  present  showing  evidence 
that  the  turbulence  model  is  operating  in  LES  mode.  Based  on  these  assessments,  the 
dimensional  time  step  was  selected  to  be  equivalent  to  a  1°  azimuthal  sweep  of  the  rotating 
hub  for  all  simulations.  Up  to  forty  subiterations  were  required  to  achieve  time  accuracy 
through  residual  reduction  by  two  orders  of  magnitude  during  each  time  step. 
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(a)  Wake  location  1 


(b)  Wake  location  2 


(c)  Wake  location  3 


(d)  Wake  location  4 


Figure  26:  Comparison  of  power  spectral  density  (PSD)  predicted  using  two  different  time 
steps  for  the  Re  =  0.29  x  106  rotating  hub  at  /i  =  0.152.  The  wake  locations  (1-4)  correspond 
to  Fig.  44.  (Note:  Data  arc  presented  using  the  final  adapted  grid.) 


6.4  Static  and  Dynamic  Force  Correlations 

The  first  simulations  undertaken  were  to  evaluate  the  grid  requirements  and  the  ability  of  the 
computational  simulations  to  predict  the  hub  forces.  With  the  requisite  number  of  normal 
cell  layers  within  the  boundary  layer,  as  discussed  in  Section  5.1.1  and  Ref.  51,  the  static 
hub  drag  is  predicted  very  accurately  by  the  computational  approach,  as  depicted  in  Fig.  27. 
The  drag  is  linear  with  the  squared  free  stream  speed.  The  maximum  error  occurs  near  a 
speed  of  zero  with  an  error  of  5%.  Most  errors  are  within  l%-2%  over  the  speed  range. 
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Figure  27:  Correlation  between  experimental  load  cell  and  CFD  data  for  the  0°  static  hub 
model 


Comparison  of  the  drag  coefficient  for  the  0°  and  45°  static  cases  and  the  rotating  case  at 
24U  rpm  results  in  drag  coefficient  values  that  are  similar.  The  maximum  drag  occurs  at  0° 
static  orientation,  while  the  minimum  occurs  at  the  45°  static  orientation.  The  drag  changes 
are  directly  related  to  the  minimization  of  the  frontal  flat  plate  area  as  the  orientation 
changes.  This  change  results  in  fewer  high  pressure  (near  stagnation)  areas,  that  arc  instead 
now  suction  areas,  as  illustrated  in  Fig.  28.  Not  unexpectedly,  the  rotational  hub  drag  and 
pressure  coefficient  behavior  falls  between  these  two  static  results.  Observing  the  pressures 
at  the  0°  azimuth  (Figs.  28  a  and  c),  overall  the  component  pressure  distributions  appear 
to  be  similar,  with  the  exception  that  the  areas  where  lower  pressure  occurs  have  higher 
magnitudes  when  the  model  is  rotating.  This  is  particularly  true  for  the  forward  and  aft 
(with  respect  to  the  incoming  flow  from  the  left)  root  struts  and  blocks.  The  additional  local 
velocity  due  to  the  rotation  (V  —  Q.R)  is  particularly  pronounced  as  the  distance  from  the 
center  of  rotation,  /?,  is  increased,  as  expected. 

Table  7:  Comparison  of  global  (total)  drag  coefficient. 


Model  Scale 


Static  hub  at  0°  orientation 

1.26 

Static  hub  at  45°  orientation 

1.19 

Rotating  hub  at  240  rpm 

1.23 
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cp:  1.5  1  -0.5  0  0.5  1 


cp  -1.5  -1  -0.5  0  0.5  1 


cp: 


-1.5  -1  -0.5  0  0.5  1 


(a)  Static  hub  at  0°  orientation  (b)  Static  hub  at  45°  orientation  (c)  Rotating  hub  at  240  rpm 


Figure  28:  Surface  pressure  coefficient  distribution  on  the  static  and  rotating  rotor  hub.  Free 
stream  flow  enters  normal  to  the  figure. 

6.5  Drag  Deconstruction 

Analysis  of  the  change  in  the  individual  component  contributions  to  drag  provides  some 
insight  into  the  cause  of  the  differences  in  the  total  drag,  and  are  directly  related  to  the 
surface  pressure  changes  in  Fig.  28.  The  drag  from  each  of  the  components  for  both  the 
static  eases  (Table  8)  and  the  rotating  cases  (Tabic  9)  arc  presented.  Using  Fig.  21  to  aid 
in  the  analysis,  the  differences  in  the  physics  is  apparent.  The  components  on  the  windward 
(forward)  side  of  the  hub  maintain  the  higher  drag,  while  the  leeward  (aft)  hub  components 
have  significant  reductions.  The  strut  for  the  45°  azimuth  has  a  higher  drag  value,  but  this 
is  offset  by  the  significantly  lower  shank  and  pitchlink  drag  values. 

The  actual  change  in  drag  coefficient  on  each  component  (Figs.  29a  and  30a)  yields 
the  overall  drag  change,  but  the  importance  of  each  change  is  more  readily  observed  by 
the  percentage  drag  change  in  Figs.  29f  and  30b.  Approximately  half  of  the  components 
contribute  a  small  relatively  constant  percentage  change  in  the  drag  (within  2%  -  4%)  no 
matter  what  the  configuration,  static  or  rotating.  For  the  components  that  arc  fiat  plates  of 
various  shapes,  the  drag  reduction  can  be  attributed  primarily  to  the  viscous  (skin  friction) 
drag  as  Reynolds  number  increases.94 

Four  components  in  particular,  the  blocks,  driveshaft,  shanks,  and  pitchlinks,  do  not  con¬ 
form  to  the  relatively  constant  trend.  These  components  are  notably  comprised  of  primarily 
bluff  bodies,  three  of  which  are  circular  cylinders.  Further  analysis  of  the  cylindrical  compo¬ 
nents  is  compared  with  experimental  results  from  Hocrncr.95  Due  to  the  installation,  most 
of  these  cylindrical  components  have  “wall-like"  ends  and  are  of  sufficient  aspect  ratio  that 
they  can  be  considered  as  essentially  two-dimensional  bodies,  as  per  Hoerner’s  experimental 
assessment.  The  critical  Reynolds  number  (Rcc)  for  cylinders  in  a  cross- flow  is  approximately 
300,000-400,000.  The  driveshaft  is  below  Rec  and  correlates  with  experimental  drag  predic¬ 
tions  of  1.18  very  well.  Pitchlink  4  most  accurately  replicates  the  cylinder  drag  prediction, 
which  is  attributed  to  its  orientation  where  it  encounters  minimal  interference  effects  and 
nearly  unperturbed  free  stream  flow.  Shank  interference  causes  a  reduction  in  drag  obtained 
for  pitchlinks  1  and  2*  pitchlink  1  is  affected  by  the  presence  of  its  shank  and  the  assembly 
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Table  8:  Drag  tabulation  for  cylindrical  components  static  azimuth  orientations. 


0° 

45° 

Component 

Rc-d 

CD 

Cd 

Driveshaft 

15,200 

1.2305 

1.612 

Pitchlink  1 

7,600 

0.9353 

0.9840 

Pitchlink  2 

7,600 

1.0338 

0.5697 

Pitchlink  3 

7,600 

0.5377 

1.1620 

Pitchlink  4 

7,600 

1.1988 

0.9614 

Shank  1 

21.000 

0.2566 

0.8151 

Shank  2 

21.000 

0.7091 

0.0469 

Shank  3 

21,000 

0.0421 

0.6929 

Shank  4 

21.000 

0.6127 

0.2639 

Table  9:  Average  drag  tabulation  for  rotating  cylindrical  components. 


Component 

Driveshaft 

Pitchlinks 

Shanks 

Rotating  at  240  rpm 

1.2172 

1.0371 

0.6243 

of  pitchlink  2  is  on  the  leeward  side  of  its  shank.  Pitchlink  3  shows  interference  effects  due 
to  its  location  farthest  downstream  with  respect  to  other  components.  Overall,  the  influence 
of  scaling  of  the  pitchlinks  is  negated  due  to  the  compensating  interference  effects,  and  this 
was  observed  at  all  static  locations.  Shanks  2  and  4,  which  should  nominally  compare  with 
Hoerner  at  drag  values  of  1  18,  encounter  interference  and/or  finite  aspect  ratio  effects  that 
result  in  reduced  drag  from  theory.  For  the  rotating  cases,  the  drag  increases  as  the  rotation 
speed  increases,  indicating  the  possible  presence  of  the  "Magnus”  effect,95  as  well  as  changes 
in  the  friction  drag  and  interference  drag. 

The  above  computations  include  the  interference  effects  on  the  components.  When  each 
component  is  removed  and  the  simulation  recomputed,  Fig.  31  indicates  that  the  pitch  links 
and  blocks  have  significant  decrease  due  to  shank  removal.  Removal  of  the  hub  assembly 
results  in  a  drag  increase  on  the  pitch  links  and  a  similar  decrease  on  the  driveshaft.  Notably, 
the  pitch  link  removal  results  in  small  change  to  surrounding  components. 
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Figure  29:  Component  drag  contributions  to  global  hub  drag  for  static  simulations.  (Note: 
Reynolds  numbers  arc  abbreviated  such  that  ‘M  implies  millions.) 
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Figure  30:  Component  drag  contributions  to  global  hub  drag  for  rotating  simulations.  (Note: 
Reynolds  numbers  are  abbreviated  such  that  ‘M'  implies  millions.) 
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Figure  31:  Model  scale  ( Re  =  0.29  M  at  /r  =  0.152)  drag  build-up  due  to  sequential  decon¬ 
struction. 


59 


Figure  32:  Comparison  of  linear  and  CFD  computed  drag  reduction  due  to  hub  components 
at  Re  =  0.29  M  and  / 1  =  0.152. 


Figure  33:  Total  drag  comparison  at  Re  =  0.29  M  and  p  =  0.152. 


Analysis  of  surface  pressure  contours  on  the  hub  may  provide  further  insight  in  to  the 
trends  in  interference  effects  (Fig.  32)  and  the  total  drag  comparison  (Fig.  33).  Out  of  the 
three  deconstruction  steps,  the  shank  removal  shows  the  biggest  departures  for  linear  drag 
reduction.  For  the  Re  =  0.29  M  condition  (Fig.  34  (a)  v.  (b)),  the  major  areas  of  drag 
reduction  arc  observed  on  the  blocks  (retreating  side),  and  the  pitch  links,  mainly  due  to 
shank  interference. 

Small  differences  arc  visually  noted  due  to  hub  removal  (sec  Fig.  34  (b)  v.  (c)).  These 
are  on  the  upper  drive  shaft,  the  hub  mount,  and  the  upper  portion  of  the  pitch  links.  Link 
removal  creates  the  largest  change  on  the  surface  pressures  resulting  on  the  drive  shaft. 
Minimal  difference  is  noted  for  the  pressures  for  the  p  =  0.152  condition  at  Re  =  0.29  M 
(Fig.  34  (c)  v.  (d)). 
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Pressure  coefficient  distributions  at  various  cross  sections  of  the  driveshaft  (sec  Fig.  35) 
are  compared  to  note  the  influence  of  interference  effects  in  Fig.  37.  The  pressure  coefficient 
is  plotted  as  a  function  of  the  central  angle  6  as  defined  in  Fig.  36.  Asymmetry  in  the 
distribution  is  expected  with  separated  flow  dominating  the  pressure  side  (6  >  0).  The 
overall  trend  for  all  conditions  with  component  removal  causes  the  asymmetry  to  decrease 
at  the  model  scale.  At  the  Re  =  0.29  M  condition  (Figs.  37  and  Figs.  38),  the  links  removed 
condition  recovers  nearly  all  the  symmetry  except  for  the  small  rotation  effects  that  alter 
the  flow  in  a  more  localized  fashion  and  are  not  affected  by  pitch  links,  the  hub  plates,  or 
the  shanks.  The  nearly  symmetric  pressure  coefficient  curves  on  the  model  scale  driveshaft 
are  not  reproduced  at  the  full  scale  simulations.96 

The  separation  point  on  the  suction  side  moves  more  upstream  at  the  extreme  cross- 
sections  (slices  one  and  eight)  for  the  full  hub  (Fig.  37  and  Fig.  38).  For  all  three  subsequent 
deconstruction  configurations,  the  separation  point  at  slice  one  is  slightly  upstream  in  com¬ 
parison  to  the  other  cross-sections.  This  pattern  is  not  observed  for  slice  eight,  on  the  other 
hand,  as  the  flow  remain  more  attached  than  at  other  cross  sections.  This  may  be  because 
slice  eight  is  affected  by  interference  from  a  static  object  (the  fuselage),  as  opposed  to  slice 
one  being  affected  primarily  by  the  rotating  hub  plates  and  mount. 


(c)  Hub  Removed  (d)  Links  Removed 


Figure  34:  Surface  pressure  contours  at  Re  =  0.29  M  and  //  =  0.152. 
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Figure  35:  Driveshaft  slice  locations  at  which  pressure  coefficient  plots  arc  investigated. 


90° 


Figure  36:  Definition  of  the  driveshaft  central  angle  6  with  respect  to  circular  components 
in  the  flow  field. 
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(a)  Full  Hub  (b)  Shanks  Removed 
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Figure  37:  Pressure  distribution  (top  four  slices)  at  Re  =  0.29  M  and  p,  =  0.152. 
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(a)  Full  Hub 


(b)  Shanks  Removed 


(c)  Hub  Removed 


(d)  Links  Removed 


Figure  38:  Pressure  distribution  (bottom  four  sliecs)  at  Re  =  0.29  M  and  /jl  =  0.152. 


The  interference  effects  between  the  components  is  further  elucidated  by  plotting  the 
pressure  coefficients  along  each  slice  for  the  deconstruction,  as  provided  in  Fig.  39.  At  the 
model  scale  Reynolds  number  (Fig.  39),  the  removal  of  the  shanks  has  a  dramatic  effect  on  the 
pressure  coefficient  at  all  slice  locations.  In  particular,  the  primary  cause  of  asymmetry  in  the 
pressure  distributions  is  clearly  caused  by  the  the  rotating  shanks.  Further  deconstruction 
results  in  minor  changes  in  the  pressure  distributions,  which  appear  to  be  very  similar  to 
infinite  cylinder  pressure  distributions.95  97 

The  specifics  of  the  pressure  distributions  arc  quantified  in  Tables  10  to  13. 
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Fun  Hub 
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(b)  Slice  2 


(c)  Slice  3 


(d)  Slice  4 


(e)  Slice  5 


(f)  Slice  6 


Figure  39:  Pressure  distributions  along  slicc6ipcations  (Fig.  35)  illustrating  component  in¬ 
terference  at  Re  —  0.29  M  and  fi  =  0.152. 


Table  10:  Full  Configuration  (Re  =  0.29  M  at  //  =  0.152) 


Slice 

Vmax 

^ Pm  ax 

^ PminT 

6b 

'"'Pmin/ 

0Cv 

rrninr 

$  sepi 

Q-Sepr 

1 

39.027 

0.731 

-4.663 

-0.709 

-80.216 

178.205 

-113.134 

66.863 

2 

15.539 

1.437 

-6.571 

-0.225 

-87.687 

81.712 

-121.576 

41.689 

3 

14.086 

1.390 

-6.610 

-0.405 

-90.000 

79.696 

-120.818 

122.041 

4 

13.491 

1.168 

-6.495 

-0.235 

-95.498 

95.556 

-126.927 

39.113 

5 

16.403 

1.059 

-6.795 

-0.124 

-101.725 

102.475 

-177.230 

60.462 

6 

14.191 

1.092 

-6.704 

-0.115 

-96.095 

83.904 

-165.807 

132.760 

7 

15.042 

0.938 

-6.737 

-0.373 

-97.888 

97.687 

-171.467 

122.847 

8 

13.036 

0.948 

-4.774 

-1.048 

-98.543 

94.378 

-177.603 

53.284 

Table  11:  Shanks  Removed  (Re  =  0.29  M  at  //  =  0.152) 


Slice 

°cv 

L Pmax 

^  Pmax 

^  Pmini 

^  Pminr 

(L 

Pmin/ 

6L 

Pmirir 

@sepi 

^ sepr 

i 

-2.783 

0.923 

-1.486 

-1.550 

-69.169 

75.598 

-118.481 

89.999 

2 

1.207 

1.018 

-1.687 

-1.610 

-71.185 

71.866 

-131.220 

85.897 

3 

0.583 

0.999 

-1.861 

-1.464 

-68.873 

74.259 

-134.957 

144.949 

4 

0.448 

0.996 

-2.297 

-1.510 

-77.113 

76.768 

-179.551 

89.999 

5 

2.768 

1.022 

-1.751 

-1.165 

-77.524 

68.720 

-133.410 

84.313 

6 

0.575 

1.021 

-1.629 

-1.486 

-74.076 

69.579 

-133.820 

83.904 

7 

-0.338 

1.010 

-1.714 

-1.522 

-73.108 

71.858 

-144.815 

89.999 

8 

2.395 

1.165 

-2.562 

-2.396 

-90.000 

89.999 

-112.271 

111.505 

Table  12:  Hub  Removed  (Re  =  0.29  M  at  //  =  0.152) 


Slice 

^mai 

^Pmai 

^ P-mini 

^  Pminr 

0Cv 

Pmin/ 

6, 

Pmvnr 

@sepi 

@sepr 

1 

-0.359 

0.979 

-0.899 

-0.679 

-65.601 

64.521 

-140.362 

75.827 

2 

1.361 

1.002 

-0.992 

-0.705 

-65.160 

65.752 

-115.475 

77.801 

3 

2.284 

1.013 

-1.188 

-0.951 

-68.450 

69.556 

-161.033 

79.792 

3 

0.179 

1.026 

-1.685 

-1.561 

-75.131 

78.309 

-179.819 

89.999 

5 

0.013 

1.039 

-1.442 

-1.180 

-72.682 

69.068 

-141.958 

83.802 

6 

0.890 

1.035 

-1.158 

-0.934 

-69.891 

69.273 

-134.683 

77.110 

7 

1.428 

1.045 

-1.576 

-1.250 

-72.693 

75.076 

-148.401 

82.859 

8 

2.291 

1.176 

-2.305 

-2.209 

-83.349 

89.999 

-177.711 

111.591 
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Table  13:  Links  Removed  (Re  =  0.29  M  at  p  —  0.152) 


Slice 

0, 

l'Pmai 

^Pmax 

r 

'Pm  mi 

^ Pmtrir 

CPmini 

K  . 

rrTHTlr 

$sepi 

@sepr 

1 

1.572 

1.001 

-1.269 

-1.243 

-68.137 

70.137 

-117.443 

80.659 

2 

-0.860 

1.009 

-1.384 

-1.502 

-68.164 

71.865 

-121.576 

85.897 

3 

0.583 

1.011 

-1.416 

-1.622 

-68.874 

74.259 

-96.936 

89.999 

4 

0.448 

0.976 

-1.712 

-1.699 

-77.113 

78.185 

-179.551 

89.999 

5 

2.768 

1.045 

-1.368 

-1.249 

-69.539 

72.440 

-155.186 

84.313 

6 

0.575 

1.040 

-1.474 

-1.139 

-71.939 

69.579 

-179.424 

83.904 

7 

1.985 

1.030 

-1.715 

-1.435 

-73.108 

71.858 

-178.014 

89.999 

8 

2.392 

1.165 

-2.456 

-2.371 

-90.000 

89.999 

-112.215 

111.505 

6.6  Interference  Drag 

Component  drag  deficits  for  the  static  models  lie  between  4%  and  54%.  The  individual 
component  analysis  does  not  account  for  these  large  deficits,  thus  interference  drag  is  clearly 
present.  The  presence,  or  lack  thereof,  of  interference  drag  can  be  readily  observed  using  sur¬ 
face  pressure  coefficients.  Figure  28  compares  the  surface  pressures  for  the  0°static,  450static, 
and  rotating  configurations,  respectively.  Pressure  contours  that  remain  consistent  along  a 
component  indicate  minimal  interference  effects  are  present.  This  permits  the  application  of 
component  drag  characteristics,  such  as  those  found  in  Hocrner95  during  design. 

The  nonlinearity  associated  with  interference  drag  for  the  rotating  hubs  is  clearly  illus¬ 
trated  in  Fig.  43.  Here,  the  rotating  cases  are  examined  in  planform  images  for  the  mid-plane 
of  the  hub  using  identical  contour  levels.  Overall,  greater  mixing  is  observed  in  the  full  scale 
flow  fields  (hence  vorticity  levels  are  lower)  due  to  increase  in  the  Reynolds  number  creating 
a  more  turbulent  wake.  The  change  in  pressure  coefficient  and  shed  wake  are  clearly  related 
when  comparing  Figs.  43  and  28.  Interference  effects  can  be  observed  on  the  left  blade  shank 
connector  (just  below  the  bolts),  in  Fig.  28  where  the  trailing  edge  pressure  suction  shifts  to 
a  more  forward  location  from  the  bottom  to  top  of  the  connector. 

The  removal  of  the  blade  shanks  to  examine  the  deconstructed  model  indicates  that  the 
blade  shank  contribution  to  the  total  drag  is  minimal.  Upon  removal  of  the  blade  shanks  and 
the  hub  plates,  the  drag  decreases  by  approximately  one-third.  Two  thirds  of  the  hub  drag 
arc  due  to  the  contributions  of  the  drive  shaft,  swash  plate  and  pitch  link  drag  increments. 
Rotating  each  of  these  deconstructed  models  confirms  that  there  is  little  variation  in  drag 
from  the  static  configuration.  The  pitch  link  drag  contribution  may  be  computed  from  the 
shift  in  measurements  from  the  prior  deconstructed  model. 

6.7  Wake  Characterization 

The  unsteady  wake  of  the  hub  was  experimentally  evaluated  through  the  use  of  particle 
image  velocimetry  (PIV),  and  was  simultaneously  or  a  priori  analyzed  with  computational 
simulations.  Experimental  velocity  deficits  in  the  wake  were  obtained  at  a  distance  one  hub 
diameter,  as  measured  from  the  hub  eenterpoint  (Fig.  40). 
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(a)  Top  view 


(b)  Isometric  View 

Figure  40:  Map  of  PIV  data  collection  plane  comprised  of  overlapping  stitches.  From  Ref.  91. 


Computational  results  were  correlated  with  experimental  data  to  ascertain  the  ability  of 
the  computational  adapted  mesh  -  obtained  a  priori  to  observing  the  experimental  results 
to  capture  the  details  of  the  wake.  Both  experimental  and  computational  methods  show 
a  contraction  of  the  momentum  deficit  in  the  hub  wake  when  the  model  is  oriented  at  45° 
azimuth,  corresponding  to  a  reduction  of  hub  drag  based  on  frontal  area.  The  experimental 
and  computational  wake  velocity  deficits  correlated  very  well,  as  illustrated  by  Fig.  41. 
The  causal  wake  behavior  that  drive  the  momentum  deficits  observed  in  Fig.  41  is  clearly 
illustrated  by  the  computational  velocity  contours  downstream  of  the  hub  (Fig.  42).  The 
largest  deficit  at  the  static  hub  centerline  is  clearly  defined  by  the  strong  velocity  deficit  just 
behind  the  hub  main  shaft.  The  two  secondary  wake  deficits  appear  behind  the  pitch  links 
to  the  left  and  right  of  the  hub  shaft.  When  the  hub  rotates  in  a  counter  clockwise  direction 
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the  primary  veloeity  deficit  is  translated  upward  and  to  the  right,  appearing  behind  the  right 
(aft  looking  forward)  blade  shank.  The  veloeity  defieit  due  to  the  main  drive  shaft  appears 
to  have  eoaleseed  with  that  of  the  right  piteh  link,  leveling  only  a  secondary  veloeity  deficit 
from  the  left  piteh  link,  whieh  has  also  translated  in  the  positive  y-direction  (to  the  right) 
of  the  flow  indueed  by  the  counter  cloekwuse  rotation. 


Figure  41:  Comparison  of  experimental  (PIV)  and  a  pnon  computational  data  for  the  tunnel 
axis  wake  veloeity  deficit  for  —  8.941  m/s. 
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(b)  Counter  clockwise  rotating  hub  at  240  rpm 


N 

-0.2 


y 

(a)  Static  hub  at  0°  azimuth  orientation 


Figure  42:  Velocity  contours  one  hub  diameter  downstream  of  the  hub  center  for  =  8.941 
m/s. 


Further  insights  into  the  character  of  the  flow  on  and  behind  the  hub  can  be  obtained 
from  the  computational  simulations.  Consider  for  example,  the  velocity  contours  presented 
in  planform  view  (looking  down  on  the  rotor  hub)  in  Fig.  43  for  the  static  45°  azimuth, 
rotating,  and  static  0°  azimuth  cases.  (The  solid  black  line  on  the  right  of  each  figure 
indicates  the  downstream  location  of  the  velocity  contours  presented  in  Fig.  42.)  The  effects 
of  rotation  yield  a  contracted  wake  span  similar  to  the  45°  static  case.  As  noted  previously 
in  the  discussion  of  the  force  measurements,  the  drag  for  the  rotating  hub  is  similar  to  the 
drag  of  the  static  45°  azimuth  orientation.  The  near-body  wake  behavior  of  these  two  eases 
are  clearly  similar,  in  contrast  to  the  near-body  shed  wake  of  the  static  0°  azimuth  case.  The 
rotation  causes  significant  vortex  shedding  from  the  blade  shanks,  similar  to  the  bluff-body 
shed  w’ake  of  the  shanks  wThen  they  are  obliquely  aligned  with  the  free  stream. 
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Figure  44:  Wake  hot-wire  sampling  locations.  From  Ref.  91. 


Figure  43:  Illustration  of  the  wake  extent  (via  vorticity  contours)  for  the  static  at  45° 
orientation  (left),  rotating  (center)  and  0°  static  orientation  (right)  (U^  =  8.941  m/s). 


6.8  Velocity  Frequency  Spectra 

Experiments  were  conducted  wherein  a  single  axis  hot-film  anemometer  measured  velocity 
fluctuations  at  a  few  selected  locations  (Fig.  44)  in  the  hub  wake.  Experiment  and  computa¬ 
tional  power  density  spectra  plots  for  the  static  and  rotating  tests  are  illustrated  in  Figs.  45 
and  46,  which  correspond  to  locations  moving  inboard  of  the  hub  wake.  The  initial  computa¬ 
tional  simulations  were  computed  before  the  experimental  data  were  collected.  Because  the 
hot-wire  anemometry  locations  were  not  known  until  after  the  experiment,  these  runs  were 
continued,  without  user  intervention,  until  enough  data  at  these  locations  were  obtained  to 
perform  the  Fast  Fourier  Transforms  (FFT)  analyses  at  the  selected  locations.  During  these 
additional  rotor  revolutions,  no  significant  changes  in  the  other  parameters  analyzed  were 
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observed. 

The  power  spectra  data  are  compared  with  the  classic  5/3  law,  which  characterizes  ho¬ 
mogeneity  of  turbulence.  A  study  in  the  shedding  frequencies  of  bluff  bodies  by  Sakamoto 
and  Arie98  provides  Strouhal  numbers  of  circular  cylinders  of  varying  aspect  ratios.  Their 
investigation  yields  predicted  shedding  frequencies  in  agreement  with  the  blade  shank  fre¬ 
quency  of  approximately  50  Hz  as  seen  from  sample  location  two,  just  behind  the  shank. 
There  is  a  noticeable  amplitude  spike  for  all  sampling  locations  of  the  rotating  test  case  cor¬ 
responding  to  the  predicted  four-per-revolution  signal  of  the  four-bladed  hub  model.  Both 
the  experimental  and  CFD  efforts  capture  this  corresponding  frequency  of  16  Hz  as  the  hub 
rotates  at  240  rpm  (4  rev/sec).  The  next  harmonic  of  this  frequency  (8  per  rev.)  is  greatly 
amplified  at  the  inboard  sampling  location,  well  within  the  wake  region. 


Figure  45:  Power  density  spectrum  of  velocity  fluctuations  in  the  hub  wake  for  the  static 
hub  at  0°  orientation  at  =  13.41  m/s. 
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Figure  46:  Power  density  spectrum  of  velocity  fluctuations  in  the  hub  wake  for  the  rotating 
hub  at  240  rpm  =  8.941  m/s. 


Figure  47:  Top  view  of  HRLES  blending  function  contours. 

Better  agreement  with  5/3  law  slope  is  noticed  for  both  the  experimental  and  compu¬ 
tational  results  at  inboard  locations.  The  turbulence  displays  an  expected  broad  range  of 
scales  in  the  bluff  body  wake  region.  Contours  of  the  IIRLES  turbulence  model’s  blending 
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function  are  given  in  Fig.  47.  Values  near  one  indicate  that  the  RANS  formulation  is  in 
effect  and  values  near  zero  indicate  regions  dominated  by  LES.  The  four  locations  show  that 
the  methodology  is  capturing  of  the  broad  range  of  turbulent  scales  with  LES. 

The  power  spectral  plots  (Figs.  45,  46,  and  48)  show  some  disagreement  in  the  magnitude 
of  power  obtained  between  experiment  and  computations.  The  general  observation  is  an  over- 
prediction  in  power  obtained  by  CFD  at  location  one,  attributed  to  the  computational  shear 
layer  prediction  coinciding  with  location  one.  One  explanation  for  the  discrepancies  are  that 
the  wake  comparisons  have  not  been  made  at  the  same  location  due  to  an  experimental 
uncertainty  (as  noted  by  the  rapid  changes  in  the  wake  extent  near  the  comparison  location 
in  Fig.  43).  Another  is  that  further  additional  adaptations  or  simulation  length  may  be 
needed  in  the  computations.  Recent  experiments  with  hot-wire  anemometry  were  analyzed, 
raising  additional  questions,  as  referenced  in  Section  7.3. 


(a)  Location  1  (b)  Location  2 


Figure  48:  Power  density  spectrum  of  velocity  fluctuations  in  the  hub  wake  for  the  rotating 
hub  at  240  rpm  at  —  22.35  m/s. 


6.9  Reynolds-averaged  Scaling  Analysis 

Under  the  final  funding  profile,  the  computational  scaling  analysis  from  wind  tunnel  model 
to  full  helicopter  was  not  possible.  However,  additional  funding  was  obtained  from  the 
Vertical  Lift  Consortium  (VLC)  to  perform  this  study.  This  analysis  is  documented  in 
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Refs.  93  and  96,  and  are  not  presented  in  this  report.  The  analysis  used  the  original 
computations  described  in  the  prior  sections  for  correlation.  The  full-scale  configuration  was 
examined  at  flight  conditions  based  on  the  UH-60A  C8534  counter"  to  provide  an  order  of 
magnitude  increase  in  the  global  Reynolds  number  compared  to  the  model  scale  experiment. 
To  illustrate  the  Reynolds  number  effects,  two  full-scale  rotating  conditions  to  provide  a 
variation  based  on  rotating  speed  and  advance  ratio  wrcrc  also  analyzed  along  with  the  static 
hub  orientations  are  discussed.  Finally,  the  inflow  effects  due  to  the  presence  of  a  fuselage 
were  modeled  using  the  Robin  fuselage. 

From  Shcnoy  ct  al.,93  the  main  conclusions  from  this  study  arc 

1.  The  boundary  layer  and  near  body  grid  should  be  designed  for  the  viscous  characteris¬ 
tics  of  the  highest  Reynolds  number  calculation.  Computations  indicate  that  the  error 
can  increase  linearly  with  the  increase  in  Reynolds  number. 

2.  The  application  of  a  mesh  adaptation  technique  that  encompasses  overset  grids  without 
regard  to  mesh  boundaries  improves  the  near  wake  and  integrated  load  predictions,  and 
permits  an  initial  grid  to  be  applied  to  a  number  of  configurations  without  the  need  to 
develop  new  grids. 

3.  Components  that  are  bluff  bodies  lead  to  nonlinear  scaling  of  the  drag  through  a 
combination  of  Reynolds  number  scaling  effects  and  changes  in  the  interference  drag. 
For  components  that  remain  clear  of  the  shed  wakes,  theoretical  and  experimental 
estimates  of  drag  for  each  individual  component  correlate  well  with  the  computational 
results. 

4.  The  determination  of  the  interference  drag  for  rotating  bodies  must  include  an  esti¬ 
mation  of  the  translational  shift  or  bias  in  the  shed  wake  due  to  the  generation  of  side 
forces  (Magnus  effect).  The  shed  wake  translates  as  a  function  of  the  local  velocity 
on  the  component  surface  with  the  free  stream  velocity  and  is  clearly  observed  in  the 
computational  evaluations. 

5.  When  sealing  the  drag  during  rotation,  velocity  scaling  based  on  the  advance  ratio 
rather  than  the  angular  velocity  of  the  rotor  appears  to  a  more  appropriate  physical 
scaling. 

7  Scissors  Configuration  Analysis 

The  hub  analysis  presented  thus  far  has  been  based  on  computations  and  experiments  per¬ 
formed  prior  to  late  2011,  when  wind  tunnel  renovations  were  undertaken.  A  final  set  of 
analyses  with  the  addition  of  a  scissors  component  (Fig.  49)  was  undertaken.  There  was  no 
PIV  data  obtained  in  the  last  experiment,  instead  hot-wire  anemometry  was  used  to  obtain 
wake  velocities  to  augment  the  load  cell  force  data.  All  computational  results  wTere  obtained 
a  priori  to  correlation  with  the  experimental  results. 
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t/oc 

8.941  m/s 

13.41  m/s 

22.35  m/s 

120  rpm 
240  rpm 

Forces,  Velocities 

Forces 

Forces,  Velocities 

Forces 

Table  14:  Matrix  of  flow  conditions  and  analysis  methods. 


(a)  without  scissors  (b)  with  scissors 

Figure  49:  Hub  configuration  with  and  without  scissors. 


The  hot-wire  anemometry  data  were  provided  in  the  form  of  velocities  of  the  form  14  as 
well  as  perturbations  of  this  velocity.2  The  quantity  14  is  the  effective  hot-wire  measured 
velocity  given  by 

14  =  \/u2  +  kvv2  +  tc2,  (25) 

where  kv  =  0.2  as  reported  by  experimental  calibration.  The  presentation  of  the  experimental 
data  in  this  manner  permitted  a  number  of  new  statistical  analyses. 

Table  14  describes  the  different  flow  conditions  analyzed  and  the  analyses  methods  per¬ 
formed  in  this  work.  Forces  are  computed  for  each  simulation  performed,  while  velocity  fields 
(mean  and  turbulent)  were  computed  and  correlated  (wherever  applicable)  for  =  8.941 
m/s  and  =  13.41  m/s  at  240  rpm. 

7.1  Hub  Drag  Loads  Analysis 

The  effect  of  the  loads  due  to  addition  of  the  scissors  is  described  in  Table  15.  The  additional 
10%  drag  increase  at  120  rpm  predicted  by  CFI)  was  also  predicted  by  the  experiment.52 
The  net  drag  increase  at  240  rpm  varies  with  the  free  stream  speed  and  does  not  correlate 

2Retricvcd  from  http://www. adl.gatech.edu/expaero/hubdrag/  on  July  2,  2013,  and  verified  through 
September  29,  2013. 


76 


Flow  Condition 

Configuration 

cD 

%  increase 

U0 o  =  13.41  m/s  at  120  rpm 

without  scissors 
with  scissors 

1.265 

1.397 

10.6% 

f/oo  =  8.941  m/s  at  240  rpm 

without  scissors 
with  scissors 

1.263 

1.463 

15.9% 

Uoo  =  13.41  m/s  at  240  rpm 

without  scissors 
with  scissors 

1.291 

1.437 

11.3% 

Ua o  —  22.35  m/s  at  240  rpm 

without  scissors 
with  scissors 

1.237 

1.312 

6.1% 

Tabic  15:  Drag  coefficient  comparison  showing  the  effect  of  scissors  at  240  rpm. 


with  the  drag  increase  at  120  rpm.  Also,  the  drag  increase  due  to  the  scissors  appears  to 
decrease  with  increasing  advance  ratio  or  free  stream  speed.  It  is  not  determinable  whether 
the  effect  of  the  drag  increase  due  to  the  scissors  is  a  function  of  the  free  speed  primarily 
since  simulations  were  not  performed  at  the  lower  rotor  speed  at  different  free  stream  speeds. 

7.2  Vibratory  Loads  Analysis 

The  unsteady  drag  and  side  force  coefficient  harmonics  for  the  scissors  hub  are  plotted  in 
Figs.  50  and  51,  respectively.  The  effect  of  the  scissors  can  be  assessed  via  these  plots  on 
the  contribution  to  these  forces.  It  is  clear  that  the  scissors  have  a  major  contribution  in 
the  two-per-rev  harmonic  for  both  forces.  For  both  forces,  the  four-per-rev  and  eight-per-rev 
harmonics  are  prevalent  and  the  four-per-rev  energy  is  generally  greater  than  that  of  the 
cight-pcr-rcv.  The  exceptions  arc  the  drag  harmonics  at  U =  8.941  m/s,  where  the  cight- 
per-rev  energy  is  stronger  for  both  configurations.  Frequencies  greater  than  four-per-rev  are 
caused  by  aerodvnamic  wake  structures  that  impinge  on  the  oncoming  geometry:  at  = 
8.941  m/s,  the  advance  ratio  is  smallest  (/i  =  1.463)  indicating  that  the  wake  influence  and 
impingement  on  the  oncoming  blade  on  the  advancing  side  is  the  greatest.  This  is  illustrated 
by  a  qualitative  comparison  of  the  wakes  of  two  simulations  in  Fig.  52.  At  the  lower  speed, 
there  arc  significant  wake  shedding  that  impinges  on  the  advancing  strut  (top)  from  the 
forward  strut  (left)  compared  to  the  wake  shedding  for  the  higher  velocity.  These  structures 
have  a  greater  contribution  on  the  advancing  blade  higher  harmonics,  where  a  significant 
drag  contribution  is  obtained.  This  explains  the  discrepancy  observed  in  Fig.  50(a). 
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Figure  50:  Comparison  of  unsteady  drag  coefficient  (Co)  harmonics  for  a  rotating  hub  with 
scissors  at  240  rpm. 
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Figure  51:  Comparison  of  unsteady  side  force  coefficient  ( Cy )  harmonics  for  a  rotating  hub 
with  scissors  at  240  rpm. 
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(a)  C/oo  -  8.941  m/s  (/a  =  1,463)  (b)  £/<*>  =  13.41  m/s  (p  =  2.195) 

Figure  52:  Q-eriterion  eontours  of  the  wake  at  Z  =  0.0  m  at  240  rpm. 


7.3  Wake  Velocity  Correlations 

Wake  velocities  were  obtained  from  the  CFD  simulations  at  X  =  ID,  X  =  2D,  and  X  3D 
downstream  of  the  hub  for  comparing  trends  and  analyzing  veloeitv  speetra.  The  spatial 
locations  of  these  probes  are  delineated  in  53.  In  addition  to  the  hub  level  at  Z  =  0.0,  data 
are  obtained  at  wake  locations  direetly  behind  the  scissors  at  Z  —  —0.204  D  so  that  the 
effect  of  the  seissors  on  the  downstream  wake  may  be  determined. 

The  a  priori  correlation  between  the  computational  and  experimental  wake  data  is  not 
as  aeeurate  when  hot-wire  anemometry  was  utilized,  when  compared  to  the  previous  a  priori 
computational  correlation  with  the  experimental  PIV  results  obtained  without  the  seissors 
(see  Section  6.7).  A  typieal  example  of  the  new  experimental  hot-wire  and  computational 
velocities  are  compared  in  Fig.  54.  The  computational  approach  for  measuring  the  velocities 
is  the  identical  to  the  approach  applied  with  the  PIV  experimental  correlations.  Unfortu¬ 
nately,  there  was  not  PIV  and  hot-wire  anemometry  experimental  data  obtained  with  both 
hot-wire  anemometry  and  PIV  at  the  same  locations  for  comparative  analysis.  Therefore,  an 
alternate  analysis  to  examine  the  discrepancies  using  the  computational  results  comparing 
the  data  where  both  sets  of  experimental  measurements  were  obtained. 

Since  the  seissor  components  are  loeated  near  the  hub  eenter  (Fig.  49),  the  predieted 
wake  extent  from  the  seissor  components  is  hypothesized  using  engineering  first  principles  to 
remain  in  the  wake  core.  The  free  stream  velocity  near  the  tunnel  walls  should  be  recovered 
for  the  same  free  stream  conditions. 

The  computational  data  for  the  2011  (PIV)  and  2013  (hot-wire  anemometry)  tests  are 
compared  in  Fig.  55.  The  PIV  data  are  presented  in  the  loeal  streamwise  velocity,  tx,  while 
the  hot-wire  data  require  a  computed  velocity,  14,  defined  in  25.  In  the  computational 
results  for  the  configuration  ineluding  the  scissors,  the  u  component  of  velocity  compares 
very  well  with  14,  indieating  that  the  contribution  of  the  local  three-dimensional  velocity 
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Figure  53:  Wake  velocity  measurement  locations  in  space.  (Red  line  indicates  traverse  of 
each  profile). 
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perturbations  remain  in  the  core  of  the  wake.  This  is  confirmed  with  a  consistent  comparison 
of  the  local  streamwise  velocity,  u,  across  the  span  extent  of  the  tunnel.  The  superposition  of 
the  shear  layers  between  the  wake  core  and  free  stream  indicate  insignificant  contributions 
from  the  addition  of  the  scissors.  The  free  stream  outside  of  the  wake  is  also  identically 
recovered  from  both  simulations,  which  is  required  if  the  free  stream  velocities  are  identical, 
as  reported  by  the  experimental  efforts.52 


Figure  54:  Comparison  of  computational  and  experimental  wake  velocity  profiles  at  one  hub 
diameter  in  the  wake. 
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Figure  55:  Comparison  of  computational  wake  velocity  profiles  at  one  hub  diameter  in  the 
wake. 


When  the  comparable  experimental  data  are  examined  (Fig.  56),  they  exhibit  a  clear 
difference  the  wake  extent  and  free  stream  recovery.  Additionally,  the  bias  toward  the  right 
of  the  wake  due  to  rotational  effects  found  in  both  PIV  and  computations  does  not  appear 
in  the  hot-wire  anemometry  measurements.  Similar  trends  in  the  local  maxima  and  minima 
are  generally  observed  for  both  methods  but  there  is  also  appears  to  be  a  positive  (upward) 
translational  bias  in  the  new  hot-wire  anemometry  data.  Therefore,  a  conclusion  from  this 
comparison  would  be  that  the  scissors  mitigate  the  wrake  drag,  wdiich  is  counterintuitive 
to  first  principles.  As  noted  previously,  there  wTere  no  PIV  data  obtained  for  comparison 
against  the  configuration  with  scissors.  When  comparing  the  twro  computational  and  the  two 
experimental  methods  on  the  same  plot  in  Fig.  57,  the  experimental  hot-wdre  results  become 
a  clear  outlier.  All  other  wake  extents  showT  good  correlation  with  each  other,  although  the 
effect  of  the  scissors  does  effect  the  velocity  profile  in  the  core  region. 

The  cause  for  the  discrepancy  in  the  hot-ware  anemometry  comparisons  can  not  be  con¬ 
firmed  given  the  available  experimental  data.52  One  hypothesis  is  that  the  tunnel  free  stream 
velocity  may  not  have  been  identical  to  the  velocity  reported.52  This  hypothesis  is  supported 
in  Fig.  56  since  the  free  stream  velocity  is  not  recovered,  and  the  expected  wrake  core  deficit 
was  less  than  hypothesized.  If  the  hot-wdre  data  are  translated  to  be  coincident  writh  the 
free  stream  velocity  levels  of  the  other  results  (Fig.  58),  a  portion  of  the  expected  results 
and  improved  correlation  is  observed.  However,  a  larger  right  bias  of  the  wake  wTould  be  ex¬ 
pected  with  a  higher  free  stream  velocity,  so  this  may  not  be  the  source  of  the  discrepancies. 
Other  potential  sources  of  error  include  a  hub  rotation  rate  different  from  the  reported  value, 
and/or  errors  associated  wdth  the  Vh  equation  constants  that  are  typically  obtained  from  cal- 
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ibrations  directly  before  and  after  the  experiment.  It  was  not  possible  to  completely  define 
the  source(s)  of  the  discrepancies  without  full  access  to  and  analysis  of  the  experimental 
data,  which  was  not  available. 


Figure  56:  Comparison  of  experimental  data  with  different  measurement  technique  and 
addition  of  scissors. 
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Figure  57:  Overall  comparison  of  both  sets  of  CFI)  and  experimental  data. 


Figure  58:  Comparison  of  computational  and  translated  experimental  wake  velocity  profiles 
at  one  hub  diameter  in  the  wake. 


84 


7.4  Wake  Sensitivity  to  Averaging  Sample  Size 

Mean  velocity  profiles  can  be  obtained  by  time  averaging  the  wake  velocities  after  a  fully 
developed  flow  field  has  been  achieved.  A  study  of  the  convergence  (and  periodicity)  of 
the  simulations  can  be  made  by  varying  the  sampling  window  from  one  revolution  to  two 
revolutions.  Figure  59  indicates  that  there  minimal  differences  in  the  profiles  when  using 
the  ultimate  or  both  the  final  two  revolutions  of  the  simulation.  The  percent  variation  of 
the  wake  velocity  deficit  integral  is  0.06%,  0.1%,  0.09%  for  the  three  downstream  locations. 
Two  revolutions  of  velocity  data  has  been  chosen  as  the  sampling  size  for  all  time-averaged 
wake  velocities. 
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(a)  x  =  ID 


(b)  x  =  2 D 


(c)  x  =  3 D 


Figure  59:  Time-averaged  wake  velocity  comparison  at  Z  =  0.0  m  illustrating  the  effect  of 
number  of  revolutions  on  the  averaging  of  the  wake  velocity. 
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7.5  Mean  Velocity  Trends 

Time- averaged  velocity  profiles  resulting  from  the  computational  predictions  at  wake  loca¬ 
tions  x  =  ID,  x  =  2D,  and  x  =  3D  for  the  different  free  stream  velocities  at  a  rotation 
speed  of  240  rpm  is  shown  in  Figs.  60  -  62.  Several  observations  can  be  made.  The  veloc¬ 
ity  deficit  generally  decreases  (the  wake  velocity  approaches  the  free  stream  recovery)  at  a 
given  spanwise  location  as  the  downstream  distance  from  the  hub  increases.  This  decrease 
is  not  linear  with  respect  to  the  distance,  rather  most  recovery  occurs  during  one  to  two  hub 
diameters  downstream.  Additionally,  the  increase  in  downstream  distance  tends  to  diffuse 
or  smooth  the  local  extrema  in  the  core.  The  wake  bias  toward  the  right  (advancing  side) 
decreases  as  the  downstream  distance  increases,  with  the  maximum  changes  before  two  hub 
diameters  in  the  wake  is  reached.  This  wake  skew  is  apparent  by  the  larger  velocity  deficits 
on  the  wake's  advancing  side. 

As  the  free  stream  velocity  increases  (higher  advance  ratios),  the  asymmetry  of  the  wake 
diminishes  in  the  wake  locations  examined.  This  implies  that  at  higher  advance  ratios  the 
wake  skew  tends  to  disappear,  which  is  expected  since  the  effect  of  rotation  is  limited  due 
to  the  dominance  of  the  free  stream  energy. 


Figure  60:  Time- averaged  wake  velocity  comparison  at  Z  0.0  m  at  different  downstream 
locations.  The  rotor  hub  conditions  arc  at  240  rpm  and  U ^  =  8. 941  m/s. 


87 


Figure  61:  Time-averaged  wake  velocity  comparison  at  Z  =  0.0  m  at  different  downstream 
locations.  The  rotor  hub  conditions  are  at  240  rpm  and  =  13.41m/,s. 


Figure  62:  Time-averaged  wake  velocity  comparison  at  Z  —  0.0  m  at  different  downstream 
locations.  The  rotor  hub  conditions  are  at  240  rpm  and  U0 0  =  22. 35m./, s. 
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7.6  Velocity  Histograms  and  Central  Moments 

In  order  to  compare  the  distribution  of  different  effective  velocities  measured  in  the  wake 
between  the  experimental  and  computational  results,  histograms  of  the  both  data  sets  at 
different  wake  locations  have  been  compared  and  assessed.  In  Figs.  63  -  68,  the  general  shape 
of  the  distributions  show'  good  correlation  in  the  velocity  axis,  although  the  velocity  transla¬ 
tion  noted  in  Section  7.3  is  evident.  This  implies  that  the  computation  accurately  captures 
the  general  range  of  experimental  velocities  writh  a  similar  distribution.  The  experimental 
data  appear  much  smoother  than  the  computation  data  as  the  number  of  computational 
samples  is  much  larger.  The  histograms  indicate  that  the  general  trend  is  followed  between 
both  data  sets.  The  experimental  data  set  exhibits  a  shift  in  the  —  y  direction  (leftward)  by 
about  Ay  =  — 0.177). 
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Figure  63:  Histogram  comparison  at  Z  =  0.0  m  and  X  =  ID.  The  rotor  hub  conditions  are 
at  240  rpm  and  =  8.941  m/s. 
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Figure  64:  Histogram  comparison  at  Z  =  0.0  m  and  X  =  2D.  The  rotor  hub  conditions  are 
at  240  rpm  and  Ux  —  8.941  m/s. 


Figure  65:  Histogram  comparison  at  Z  =  0.0  m  and  X  =  3D.  The  rotor  hub  conditions  are 
at  240  rpm  and  U ^  8. 941m/, s. 
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Figure  66:  Histogram  comparison  at  7j  =  0.0  m  and  X  =  ID.  The  rotor  hub  conditions  are 
at  240  rpm  and  13.41m/, s. 
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Figure  67:  Histogram  comparison  at  Z  =  0.0  m  and  X  =  2D.  The  rotor  hub  conditions  are 
at  240  rpm  and  U0 0  =  13.41m/s. 
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Figure  68:  Histogram  comparison  at  Z  =  0.0  m  and  X  =  3D.  The  rotor  hub  conditions  are 
at  240  rpm  and  Ux  =  13.41m/.s. 


Central  moment  statistics  provide  additional  details  about  the  distributions.  The  stan¬ 
dard  deviation,  which  arises  from  the  second  central  moment,  shows  dispersion  of  the  data 
away  from  the  mean.  The  skewness,  the  third  central  moment,  enables  understanding  of  bias 
of  the  data  above  or  below  the  mean.  Finally,  the  kurtosis  or  distribution  slope,  the  fourth 
central  moment,  denotes  the  probability  of  extreme  events.  High  kurtosis  implies  a  higher 
probability  of  outliers  and  indicates  a  distribution  with  longer  tails.  A  Gaussian  distribution 
will  have  zero  skewness  and  a  kurtosis  value  of  three. 

Figures  69  -  74  plot  the  experimental  and  computation  distributions  of  standard  devi¬ 
ation,  skewness,  and  kurtosis.  The  values  of  standard  deviation,  skewness,  and  kurtosis 
correlate  well  between  computational  and  experiment,  especially  in  the  middle  of  the  wake 
(—0.4  <  y/D  <  0.4).  Again,  the  experimental  statistics  are  smoother  than  that  of  computa¬ 
tions,  similar  to  the  histograms.  The  skewness  and  kurtosis  plots  of  both  data  sets  indicate 
that  the  inid-wake  velocities  have  distributions  close  to  Gaussian  distribution.  The  skewness 
shows  that  the  velocities  are  slightly  biased  lower  than  the  mean  velocities.  In  common 
with  the  prior  data  analyses,  there  is  a  clear  shift  in  the  wake  extent  of  the  experimental 
data.  From  a  statistical  point  of  view,  the  two  data  sets  give  reasonable  correlation,  with 
the  only  major  differences  being  the  smoothness  of  the  data  and  the  bias  of  the  wake  in  the 
experimental  data. 


96 


Figure  69:  Comparison  of  statistical  central  moments  from  adapted  grid  at  Z  =  0.0  m  and 
X  =  ID.  The  rotor  hub  conditions  are  at  240  rpm  and  —  8.941  m/s. 
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(a)  Standard  Deviation 


Figure  70:  Comparison  of  statistical  central  moments  from  adapted  grid  at  Z  =  0.0  m  and 
X  =  2D.  The  rotor  hub  conditions  are  at  240  rpm  and  U0 0  =  8.941  m/s. 
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Figure  71:  Comparison  of  statistical  central  moments  from  adapted  grid  at  Z  =  0.0  m  and 
X  =  3D.  The  rotor  hub  conditions  are  at  240  rpm  and  U0 0  =  8.941m/s. 
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Figure  72:  Comparison  of  statistical  central  moments  from  adapted  grid  at  Z  =  0.0  m  and 
X  =  ID.  The  rotor  hub  conditions  arc  at  240  rpm  and  U ^  =  13.41m/s. 
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Figure  73:  Comparison  of  statistical  central  moments  from  adapted  grid  at  Z  =  0.0  m  and 
X  =  2D.  The  rotor  hub  conditions  arc  at  240  rpm  and  =  13.41m/s. 


101 


(a)  Standard  Deviation 


Figure  74:  Comparison  of  statistical  central  moments  from  adapted  grid  at  Z  =  0.0  m  and 
X  =  3D.  The  rotor  hub  conditions  are  at  240  rpm  and  =  13.41m/s. 
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7.7  Velocity  Spectral  Analysis 

Wake  velocity  spectra  from  experiment  and  computation  arc  plotted  and  compared  in 
Figs.  75  -  80  at  the  hub  level  (Z  =  0.0  m).  The  contours  indicate  the  power  spectral 
density  (in  logarithmic  scale).  The  frequency  axis  is  in  the  vertical  direction,  marked  by 
the  per-rev  frequency  so  that  wake  structures  may  be  traced  to  their  originating  geometry. 
Noting  that  the  contour  ranges  are  the  same  for  both  computation  and  experiment,  a  gen¬ 
eral  trend  is  observed  that  the  computational  predictions  for  the  power  spectra  content  are 
greater  than  their  experimental  counterparts.  Also,  since  the  computational  data  arc  aver¬ 
aged  for  many  more  data  points  over  the  span  and  over  two  revolutions,  the  data  appear  to 
be  noisy.  However,  clear  correlations  are  observed  particularly  in  the  data  at  the  four-per-rev 
and  two-per-rev  levels. 

The  four-per-rev  levels  structures  are  generally  present  over  the  range  the  of  the  wake, 
but  they  are  much  stronger  in  the  retreating  portion  of  the  wake  (left  side  of  the  plots).  These 
structures  emanate  primarily  from  the  blade  shanks  and  pitch  links.  These  structures  also 
persist  to  the  X  =  3D  downstream  location,  the  farthest  point  measured  by  experiment.  The 
eight-per-rev  and  sixteen-per-rev  content  are  much  stronger  for  =  13.41  m/s,  but  they 
diminish  more  rapidly  at  the  downstream  locations  compared  to  the  four-per-rev  structures. 
Most  of  these  structures  are  strongest  in  the  retreating  portion  of  the  wake,  where  the 
spanwise  mixing  of  the  flow  is  much  less  compared  to  the  advancing  rotor  wake  due  to  the 
skewed  wake  shape.  A  considerable  two-per-rev  shedding  is  also  observed  in  both  sets  of 
data  and  is  intermittently  prominent  across  the  wake. 

Additionally,  traces  of  six-per-rev  shedding  are  observed.  This  corresponds  to  a  frequency 
of  24  Hz,  which  can  be  attributed  to  the  vortex  shedding  off  the  blade  shanks  (D  =  0.0349 
m)ee.  Sakamoto  and  Arie98  have  observed  that  for  cylinders  of  aspect  ratio  2.5,  which 
corresponds  to  the  aspect  ratio  of  the  blade  shanks,  the  Strouhal  number  (St  —  ^)  can  vary 
between  0.11  and  0.14,  depending  on  the  portion  of  the  cylinder  submerged  in  the  boundary 
layer.  Using  the  tip  speed  Vtip ,  as  the  characteristic  velocity  of  shedding,  a  Strouhal  number 
of  0.135  is  obtained  for  this  configuration.  Therefore,  the  six-per-rev  structures  here  can  be 
traced  to  the  rotating  blade  root.  This  is  similar  to  the  findings  reported  by  Reich  et  al.,21 
where  a  six-per-rev  wake  structure  was  attributed  to  a  Strouhal  shedding  due  to  the  hub 
arms  of  a  four-bladed  rotor  hub. 
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Figure  75:  Wake  PSD  comparison  at  Z  =  ().()  m  and  X  =  ID.  The  rotor  hub  conditions  are 
at  240  rpm  and  U ^  =  8.941m/s. 
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Figure  76:  Wake  PSD  comparison  at  7  —  0.0  m  and  X  =  2D.  The  rotor  hub  conditions  are 
at  240  rpm  and  Ux  =  8.941m/.s\ 
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Figure  77:  Wake  PSD  comparison  at  Z  —  0.0  m  and  X  =  3D.  The  rotor  hub  conditions  are 
at  240  rpm  and  =  8.941ra/s. 
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Figure  78:  Wake  PSD  comparison  at  Z  =  0.0  m  and  X  =  ID.  The  rotor  hub  conditions  are 
at  240  rpm  and  U ^  =  13.41m/s. 
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Figure  79:  Wake  PSD  comparison  at  Z  =  0.0  m  and  X  =  2D.  The  rotor  hub  conditions  are 
at  240  rpm  and  13.41m/s. 
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Figure  80:  Wake  PSD  eomparison  at  Z  =  0.0  m  and  X  =  3 D.  The  rotor  hub  conditions  are 
at  240  rpm  and  =  13.41m/s. 


Figures  81-  83  plot  the  computational  velocity  spectra  in  the  wake  at  Z  =  — 0.204D, 
illustrating  the  wake  character  at  behind  the  scissors  geometry.  Unlike  the  wake  at  the  hub 
level  (Z  =  0.0),  the  wake  extent  here  is  narrower,  but  the  rightward  bias  is  still  observed, 
which  is  expected.  It  is  clear  that  in  addition  to  the  four-per-rev  structures,  there  arc 
also  two-per-rev  structures  that  are  not  as  prominently  at  the  hub  level.  The  two-per-rev 
structure  persists  strongly  up  until  the  X  =  3D  location  whereas  the  four-per-rev  features 
tend  to  diminish. 
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Figure  81:  Wake  PSD  comparison  at  Z  =  —0.204 D  and  X  =  ID  at  240  rpm. 
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Figure  82:  Wake  PSD  comparison  at  Z  =  —  0.204 D  and  A"  =  2D  at  240  rpm. 
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Figure  83:  Wake  PSL)  comparison  at  Z  —  — 0.204D  and  X  =  3 D  at  240  rpm. 


7.8  Long-age  Wake  Analysis 

Figure  84  shows  the  capability  of  the  current  computational  methodology.  That  is,  a  compu¬ 
tational  analysis  that  combines  an  adaptive  grid  capability  with  LES  wake  modeling  enables 
the  capture,  in  a  computationally  efficient  manner,  of  the  velocit}r  spectra  much  farther  in 
the  wake,  here  at  five  diameters  dowmstream.  Experimental  data  arc  not  available  at  loca¬ 
tions  beyond  X  =  3 D  due  to  the  limitations  of  the  test  section  dimensions.  The  flattening 
of  the  spectra  at  high  frequencies  indicates  that  further  simulation  time  is  required,  but 
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the  predictions  are  parallel  to  the  5/3  law  at  lower  frequencies  node  that  the  simulation  is 
capturing  those  turbulent  scales.  In  addition,  the  harmonics  associated  with  the  different 
components  are  also  still  visible.  For  the  Y  =  —D/2  location,  the  four-per,  six-per,  eight- 
per,  and  twelve-per  revolution  are  captured,  similar  to  that  of  ,21  There  were  no  two-per-rev 
structures  as  this  configuration  studied  did  not  include  the  scissors. 
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Figure  84:  Wake  spectra  comparison  at  Z  =  0.0  m  and  X  =  5D  at  240  rpm. 
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8  Summary  of  Findings 

Using  a  combination  of  experimental  and  computational  methods,  the  characterization  of 
the  integrated  loads  and  complex  wake  field  of  a  scaled  helicopter  main  rotor  hub  has  been 
obtained.  From  this  investigation,  it  can  be  concluded  that 

1.  Large  regions  of  separated  flow  contribute  to  considerable  drag  from  the  complex  hub. 
The  azimuthal  variation  of  the  static  hub  model  shows  variation  in  drag  corresponding 
to  the  projected  frontal  area  of  the  hub.  For  the  four-bladcd  hub  configuration,  the  45° 
azimuthal  orientation  yields  similar  drag  measurements  to  the  mean  values  obtained 
for  the  rotating  hub. 

2.  Initial  deconstruction  of  the  hub  configuration  has  provided  insight  into  the  drag  break¬ 
down  and  interference  effects  of  geometric  components.  At  the  Reynolds  numbers  and 
free  stream  velocities  examined,  the  hub  plates  and  blade  shanks  contribute  approxi¬ 
mately  1/3  of  the  hub  drag,  while  the  remaining  contribution  is  due  to  flow  separation 
about  the  drive  shaft,  pitch  links  and  swashplate  in  that  order.  Rotation  effects  on  the 
deconstructed  model  arc  observed  to  be  minimal  at  the  flow  conditions  evaluated. 

3.  Computational  prediction  at  the  fidelity  level  of  the  FUN3U  Navier-Stokes  solver,  has 
been  shown  to  be  an  invaluable  tool  to  augment  experimental  analyses  for  this  complex 
hub  configuration.  Strong  correlations  between  CFD  and  both  experiment  load  and 
PIV  data  permit  the  utilization  of  the  surface  and  flow  field  characteristics  predicted 
by  CFD  to  further  explain  and  clarify  the  causal  physics  due  to  the  complex  hub 
geometry. 

4.  New  CFD  capabilities  developed  to  permit  unstructured  overset  anisotropic  feature- 
based  adaptation  across  both  background  and  near-body  grids  are  essential  to  capture 
the  correct  physics  of  complex  configurations  where  significant  wake  interactions  occur 
in  both  the  near-  and  far-field  grids. 

5.  The  hub  characteristics  were  identified  with  computational  methods  for  both  full  as¬ 
sembly  and  deconstructions  of  the  hub.  From  these,  identification  of  where  theoretical 
approximations  can  be  applied  were  identified.  Interference  effects  were  also  quantified, 
and  the  sources  identified,  including  the  influence  of  a  fuselage. 

6.  Strouhal  shedding  and  wake  interactions  were  computed  from  the  rich  data  provided 
from  the  computational  simulations.  Important  effects  of  scaling  for  both  static  and 
rotating  hubs  were  quantified;  it  is  clear  that  there  are  significant  differences  in  the 
behavior  of  the  hub  and  the  component  interference  between  model  and  full  helicopter 
scales.  Turbulent  spectra  confirmed  the  ability  to  maintain  complex  wake  behavior 
over  these  long  periods. 

7.  Long-age  wake  data  were  captured  and  correlated  with  experimental  data  from  other 
sources  to  show  that  prediction  and  analysis  of  empennage-based  aeroelastic  and  un- 
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steady  aerodynamic  phenomena,  such  as  tail  buffet  and  “wag”  are  within  reach  of 
current  computational  resources. 

Computational  methods  have  now  reached  a  maturity  that,  if  care  is  taken  in  the  grid 
generation  and  turbulence  modeling,  it  may  be  more  cost-effective  and  accurate  to  design 
and  analyze  complex  hubs  directly  with  computational  methods  rather  than  reliance  on 
model  scale  experiments.  Advanced  turbulence  models  that  use  detached  or  large  eddy 
simulations  and  grid  adaptation  are  recommended  to  improve  the  performance  quantities 
and  the  unsteady  wake  characteristics. 


116 


9  Management 

9.1  Research  Leveraging 

This  work  leveraged  the  funding  to  study  the  scaling  effects  of  hubs  from  the  Vertical  Lift 
Consortium.96  Mr.  Rajiv  Shcnoy  also  spent  two  summers  as  a  NASA  Langley  Research 
Summer  Scholar  (LARSS)  while  working  on  portions  of  the  grid  adaptation.  Finally,  this 
research  has  leveraged  the  turbulence  modeling  from  the  DARPA  Quiet  Helicopter  program 
and  Task  3  of  the  Georgia  Tech  Vertical  Lift  Research  Center  of  Excellence  (VLRCOE). 

9.2  Technology  Transfer 

The  modifications  made  to  the  FUN3D  code  have  been  provided  to  the  NASA-Langlcv 
FUN3D  Development  group,  and  have  or  will  be  released  in  production  versions  of  the  code. 

Several  projects  are  currently  using  the  overset  grid  adaptation  developed  in  this  portion 
of  the  project.  They  include  Task  10  of  the  Georgia  Tech  Vertical  Lift  Research  Center  of 
Excellence  (VLRCOE),  “Aerodynamic  and  Dynamic  Interaction  of  Bluff  Bodies'’  and  a  new 
Army  Research  Office  project  on  Dynamic  Stall,  both  of  which  have  as  the  PD  (co-PI  with 
one  other  faculty  member  on  each)  Prof.  Smith. 

9.3  Awards 

•  Rajiv  Shcnoy,  Vertical  Flight  Foundation,  2013  PhD  Scholarship  Winner 

•  Rajiv  Shenoy,  Vertical  Flight  Foundation,  2011  MS  Scholarship  Winner 

•  Rajiv  Shenoy,  “Scaling  Hub  Drag  from  Model  to  Full  Scale,”  MSAE  Awarded  December 
2011 

9.4  Publications 

A  number  of  conference  and  journal  papers  have  resulted  from  funding  (or  partial  funding) 
from  this  effort.  These  papers  include: 

•  Shenoy,  R.,  Smith,  M.  J.,  and  Park,  M.  A.,  “Unstructured  Overset  Mesh  Adaptation 
with  Turbulence  Modeling  for  Unsteady  Aerodynamic  Interactions,”  AIAA  Journal  of 
Aircraft,  in  press,  2013. 

•  Raghav,  V.,  Shenoy,  R.,  Smith,  M.  L,  and  Komerath,  N.  M.,  “Investigation  of  Drag 
and  Wake  Turbulence  of  a  Rotor  Hub,”  Aerospace  Science  and  Technology,  Vol.  28, 
(1),  doi:  10. 1016/j.ast.2012. 10.012,  2013,  pp.  164175. 

•  Shcnoy,  R.  “An  Adaptive  Mesh  Refinement  Strategy  for  Static  and  Dynamic  Overset 
Grids,”  Overset  Grid  Symposium  Student  Poster  Competition,  October,  2012. 
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•  Ortega,  F.,  Shenoy,  R.,  Raghav,  V.,  Smith,  M.  .T.,  and  Komerath,  N.,  ‘‘Exploration  of 
the  Physics  of  Hub  Drag,”  Paper  AIAA-2012-1070,  AIAA  Aerospace  Sciences  Meeting 
and  Exhibit,  Knoxville,  TN,  .January  912,  2012. 

•  Shenoy,  R..  Raghav,  V.,  Ortega,  F.,  Smith,  M.  J.,  and  Komerath,  N.,  “Deconstructing 
Hub  Drag,”  Paper  AIAA-2011-3821,  AIAA  29th  Applied  Aerodynamics  Conference, 
June  2730,  2011. 

•  Shenoy,  R.  and  Smith,  M.  J.,' Unstructured  Overset  Adaptive  Mesh  Refinement  for 
Rotorcraft  Aerodynamic  Interactions,”  Proceedings  of  the  67th  Annual  Forum  of  the 
American  Helicopter  Society,  May  25,  2011. 

The  additional  funding  obtained  from  the  Vertical  Lift  Consortium  (VLC)  for  the  scaling 
computations  has  resulted  in  the  following  papers,  which  also  utilized  results  from  this  effort  : 

•  Shenoy,  R.,  Holmes,  M.,  Smith,  M.  J.,  and  Komerath,  N.,  “Scaling  Evaluations  on  the 
Drag  of  a  Hub  System,”  Journal  of  the  American  Helicopter  Society ,  Vol.  58,  (3),  July 
2013,  pp.  113. 

•  Shenoy,  R.,  Smith,  M.  J.,  and  Komerath,  N.,  “Computational  Investigation  of  Hub 
Drag  Deconstruction  from  Model  to  Full  Scale,”  Proceedings  of  the  37th  European 
Rotorcraft  Forum,  Gallarata,  Italy,  September  1215,  2011. 

In  addition,  Mr.  Rajiv  Shenoy  is  planning  to  defend  his  PhD  dissertation  in  December, 
2013,  which  wall  include  some  of  the  computations  reported  here.  It  will  include  additional 
computations  using  the  configurations  explored  in  this  effort.  This  report  wall  be  available 
through  the  Georgia  Tech  Smart  Tech  system:  http://wavw’.smartcch.gatcch.cdu. 

9.5  Data  Availability 

Animations  that  show  in  additional  detail  the  physics  of  the  wrakcs  of  the  hub  from  compu¬ 
tational  simulations  are  available  at  http://wwwT.msmith. gatech.edu/research/hubdrag. 

Larger  data  files,  including  CAD  geometry,  computational  grids,  and  computational  wake 
data  arc  planned  to  be  incorporated  into  the  “big  data”  resources  being  developed  at  the 
Georgia  Institute  of  Technology,  of  wdiich  Prof.  Marilyn  Smith  (the  PI  of  this  task)  is  a 
participant.  Links  to  this  web  resource  will  be  provided  on  the  afore  provided  website  or  can 
be  obtained  via  email  request  to  the  PI  of  this  task.  Prof.  Marilyn  Smith. 

9.6  Training 

This  effort  provided  additional  unfunded  training  for  Marlin  Holmes  (BSAE,  May  2013)  and 
Philip  Cross  in  the  area  of  CFD  data  reduction. 


118 


10  Acknowledgments 

Computational  support  was  provided  through  the  DoD  High  Performance  Computing  Cen¬ 
ters  at  ERDC  through  an  HPC  grant  from  the  US  Navy.  Any  opinions,  findings,  and 
conclusions  or  recommendations  expressed  in  this  material  are  those  of  the  author(s)  and 
do  not  necessarily  reflect  the  views  of  the  Department  of  the  Navy  or  the  Office  of  Naval 
Research. 

The  authors  would  like  to  acknowledge  and  thank  the  NASA  FUN3D  development  team, 
in  particular  Dr.  Michael  Park  and  Dr.  Eric  Nielsen,  who  pioneered  the  grid  adaptation 
efforts  within  FUN3D. 

The  authors  would  like  to  thank  Marlin  Holmes  and  Philip  Cross  for  their  help  in  the 
data  reduction  of  the  computational  analvsis  for  both  this  and  the  Vertical  Lift  Consortium 
efforts. 

The  authors  would  also  like  to  acknowledge  Profs.  Jerry  Seitzman  and  Ari  Glezer  of  the 
Georgia  Institute  of  Technology  who  provided  insight  and  help  in  analyzing  the  experimental 
data,  especially  the  hot  wire  anemometry  experimental  data  from  the  2013  tests. 


119 


11  References 


1]  Williams,  R.  and  Montana,  P.,  “A  Comprehensive  Plan  for  Helicopter  Drag  Reduction,” 
Technical  report,  Naval  Ship  Research  and  Development  Center,  Maryland,  1975. 

[2]  Keys,  C.,  Wiesner,  R.,  et  al.,  “Guidelines  for  Reducing  Helicopter  Parasite  Drag,” 
Journal  of  the  American  Helicopter  Society ,  Vol.  20,  1975,  pp.  31. 

[31  Hoffman,  J.,  “The  Relationship  Between  Rotorcraft  Drag  and  Stability  and  Control,” 
Proceedings  of  the  31st  Annual  National  Forum  of  the  American  Helicopter  Society, 
May  1975. 

[4]  Kerr,  A.,  “Effect  of  Helicopter  Drag  Reduction  on  Rotor  Dynamic  Loads  and  Blade 
Life,”  Proceedings  of  the  American  Helicopter  Society  Symposium  on  Helicopter  Aero¬ 
dynamic  Efficiency,  March,  1975. 

[5]  Wake,  B.,  Hagen,  E.,  Ochs,  S.,  and  Matalanis,  C.,  Assessment  of  Helicopter  Hub  Drag 
Prediction  with  an  Unstructured  Flow  Solver,”  65th  Annual  Forum  of  the  American 
Helicopter  Society,  May  2009. 

[6]  Saltzman,  E.  and  Ayers,  T.,  “A  Review  of  Flight-to-Wind  Tunnel  Drag  Correlation,' 
Paper  AIAA-1981-2475,  1981. 

[7]  Bushnell,  D.,  “Sealing:  Wind  Tunnel  to  Flight,”  Ann.  Rev.  of  Fluid  Mechanics , 
Vol.  38,  (1),  2006,  pp.  111-128. 

[8]  Vassberg,  J.  and  ct  al.,  “Summary  of  the  Fourth  AIAA  CFD  Drag  Prediction  Work¬ 
shop.”  Paper  AIAA-201 0-4547,  2010. 

[9]  Ueno,  M.,  Akatsuka,  J.,  and  Hidaka,  A.,  “Drag  Decomposition  Analysis  of  CFD  Data 
of  the  DLR.-F6  Model,”  Paper  AIAA-2008-6903,  August  2008. 

[10]  McCroskey,  J.,  “Computations  of  Unsteady  Separating  Flows  over  an  Oscillating  Air¬ 
foil,”  AIAA  . Journal ,  Vol.  35,  (7),  1997,  pp.  1235-1238. 

[11’  Keys,  C.  and  Rosenstein,  II. ,  “Summary  of  Rotor  Hub  Drag  Data,”  Technical  Report 
NASA  CR.-152080,  Ames  Research  Center,  1978. 

[12]  Potsdam.  M.  and  Le  Pape,  A.,  “CFD  Investigations  on  a  NACA0036  Airfoil  with  Active 
Flow  Control,”  Paper  AIAA-2008-3869,  4th  Flow  Control  Conference,  Seattle,  Wash¬ 
ington,  June  23-26  2008. 

[13]  Gregory,  .1.,  Porter,  C.,  and  McLaughlin.  T.,  “Circular  Cylinder  Wake  Control  using 
Spatially  Distributed  Plasma  Forcing,”  38  th  AIAA  Fluid  Dynamics  Conference  and 
Exhibit,  39  th  AIAA  Plasmadynamies  and  Lasers  Conference,  40  th  AIAA  Thermo¬ 
physics  Conference,  5  th  AIAA  Theoretical  Fluid  Mechanics  Conference,  4  th  AIAA 
Flow  Control  Conference,  and  the  26  th  AIAA  Aerodynamic  Measurement  Technology 
and  Ground  Testing  Conference,  2008. 


120 


14]  Roshko,  A.,  “Perspectives  on  bluff  body  aerodynamics,”  Journal  of  Wind  Engineering 
and  Industrial  Aerodynamics,  Vol.  49,  (1-3),  1993,  pp.  79-100. 

15]  Sheehy,  T.  and  Clark,  D.,  “A  Method  for  Predicting  Helicopter  Hub  Drag,”  Technical 
report,  United  Technologies  Corp.,  Stratford  CT,  Sikorsky  Aircraft  Div.,  1976. 

[16]  Sheehy,  T.  et  al .,  “A  general  review  of  helicopter  rotor  hub  drag  data,”  Journal  of  the 
American  Helicopter  Society ,  Vol.  22,  1977,  pp.  2. 

[17]  Ochs,  S.  S.,  Matalanis,  C.  G.,  Wake,  B.  E.,  and  Egolf,  T.  A.,  “Evaluation  of  Helios 
CFD  Toolset  for  Faired  Rotor-Hub  Drag  Prediction,”  American  Helicopter  Society  67th 
Annual  Forum  Proceedings,  Virginia  Beach,  Virginia,  May  3-5  2011. 

[18]  Bridgeman,  J.  O.  and  Lancaster,  G.  T.,  “Predicting  Hub  Drag  on  Realistic  Geometries,” 
American  Helicopter  Society  Aeromechanics  Specialists  Conference,  San  Francisco,  Cal¬ 
ifornia,  January  20-22  2010. 

[19]  Bridgeman,  J.  O.  and  Lancaster,  G.  T.,  “Physics-Based  Analysis  Methodology  for 
Hub  Drag  Prediction,”  American  Helicopter  Society  66th  Annual  Forum  Proceedings, 
Phoenix,  Arizona,  May  11-13  2010. 

^20  Hill,  M.  J.  and  Louis,  M.  E.,  “Rotating  Hub  Drag  Prediction  Methodology,”  American 
Helicopter  Society  Specialists  Conference  on  Future  Vertical  Lift  Aircraft  Design,  San 
Francisco,  California,  January  18-20  2012. 

[21]  Reich,  D.  and  Elbing.  B.  and  Berezin.  C.  and  Schmitz.  S.,  “Water  Tunnel  Flow  Diag¬ 
nostics  of  Wake  Structures  Downstream  of  a  Model  Helicopter  Rotor  Hub,”  American 
Helicopter  Society  69th  Annual  Forum  Proceedings,  Phoenix,  Arizona,  May  2013. 

[22]  Ruffin,  S.,  O'Brien,  D.,  Smith.  M.,  Hariharan,  N.,  Lee,  J.,  and  Sankar,  L.,  “Compar¬ 
ison  of  Rotor- Airframe  Interaction  Utilizing  Overset  Unstructured  Grid  Techniques,” 
Proceedings  of  the  42nd  AIAA  Aerospace  Sciences  Meeting  and  Exhibit,  Jan  2004. 

[23]  O  Brien,  D.,  and  Smith,  M.,  “Understanding  the  Physical  Implications  of  Approximate 
Rotor  Methods  Using  an  Unstructured  CFD  Method,”  Proceedings  of  the  31st  Annual 
European  Rotorcraft  Forum,  Sep  2005. 

[24]  O'Brien,  D.  M.,  Analysis  of  Computational  Modeling  Techniques  for  Complete  Rotor- 
craft  Configurations.  Ph.D.  thesis,  Georgia  Institute  of  Technology,  2006. 

[25;  Potsdam,  M.,  Smith,  M.,  and  Renaud.  T.,  “Unsteady  Computations  of  Rotor-Fuselage 
Interactions,”  Proceedings  of  the  35th  Annual  European  Rotorcraft  Forum,  Sep  2009. 

[26]  Smith,  M.J.,  Shenoy,  R.,  Kenyon,  A.R.,  and  Brown,  R.E.,  “Vorticity  Transport  and 
Unstructured  RANS  Investigation  of  Rotor-Fuselage  Interactions,”  Proceedings  of  the 
35th  Annual  European  Rotorcraft  Forum,  Sep  2009. 


121 


[27]  Komerath,  N.,  Smith,  M.  J.,  Tung,  C.,  “A  Review  of  Rotor  Wake  Physics  and  Model¬ 
ing,”  Journal  of  the  American  Helicopter  Society ,  Vol.  to  appear,  2011. 

[28]  Kang,  H.J.  and  Kwon,  O.J.,  “Effect  of  Wake  Adaptation  on  Rotor  Hover  Simulations 
Using  Unstructured  Meshes,”  Journal  of  Aircraft,  Vol.  38,  (5),  2001. 

[29]  Kang,  H.J.  and  Kwon.  O.J.,  “Unstructured  Mesh  Navier-Stokes  Calculations  of  the 
Flow  Field  of  a  Helicopter  Rotor  in  Hover.”  Journal  of  the  American  Helicopter  Society , 
Vol.  47,  (2),  2002. 

[30]  Dindar,  M.,  Shepherd,  M.,  Flaherty,  J.,  and  Jansen,  K.,  “Adaptive  CFD  analysis  for 
rotorcraft  aerodynamics,”  Computer  Methods  in  Applied  Mechanics  and  Engineering, 
Vol.  189,  (4),  2000. 

[31]  Potsdam,  M.  and  Mavriplis,  D.,  “Unstructured  Mesh  CFD  Aerodynamic  Analysis  of 
the  NREL  Phase  VI  Rotor,”  47th  AIAA  Aerospace  Sciences  Meeting,  Jan  2009. 

[32]  Park,  Y.M.  and  Kwon,  O.J.,  “Simulation  of  Unsteady  Rotor-Fuselage  Interactions  Using 
Unstructured  Adaptive  Sliding  Meshes,”  Journal  of  the  American  Helicopter  Society , 
Vol.  49,  (4),  2006. 

[33]  Cavallo,  P.A.,  Sinha,  N.,  and  Feldman,  G.M.,  “Parallel  Unstructured  Mesh  Adaptation 
Method  for  Moving  Body  Applications,”  AIAA  Journal  Vol.  43,  (9),  2005. 

[34]  Meakin,  R.L..  “On  adaptive  refinement  and  overset  structured  grids,”  Proceedings  of 
the  13rd  AIAA  Computational  Fluid  Dynamics  Conference,  June  1997. 

[35]  Ilenshaw,  W.D.,  and  Sclrwendeman,  D.W.,  “Parallel  Computation  of  Three-dimensional 
Flows  Using  Overlapping  Grids  with  Adaptive  Mesh  Refinement,”  Journal  of  Compu¬ 
tational  Physics ,  Vol.  227,  (16),  2008. 

136  Kannan,  R.,  and  Wang,  Z.J.,  “Overset  Adaptive  Cartesian/Prism  Grid  Method  for 
Stationary  and  Moving-Boundary  Flow  Problems,”  AIAA  Journal,  Vol.  45,  (7),  2007. 

[37]  Canonne,  E.,  Benoit,  C.,  and  Jeanfaivre,  G.,  “Cylindrical  mesh  adaptation  for  isolated 
rotors  in  hover,”  Aerospace  Science  and  Technology ,  Vol.  8,  (1),  2004. 

[38]  Duque,  E.P.N.,  Biswas,  R.  and  Strawn,  R.C.,  “A  Solution  Adaptive  Struc¬ 
tured/Unstructured  Overset  Grid  Flow  Solver  with  Applications  to  Helicopter  Rotor 
Flows,”  Proceedings  of  the  13rd  AIAA  Applied  Aerodynamics  Conference,  June  1995. 

[39]  Park,  M.A.  and  Darmofal,  D.L.,  “Parallel  Anisotropic  Tetrahedral  Adaptation,”  Pro¬ 
ceedings  of  the  46th  AIAA  Aerospace  Sciences  Meeting  and  Exhibit,  Jan  2008. 

[40]  Lee-Rausch,  E.M.,  Park,  M.A.,  Jones,  W.T.,  Hammond,  D.P,  and  Nielsen.  E.J.,  “Appli¬ 
cation  of  Parallel  Adjoint-Based  Error  Estimation  and  Anisotropic  Grid  Adaptation  for 
Three-Dimensional  Aerospace  Configurations,”  Proceedings  of  the  23rd  AIAA  Applied 
Aerodynamics  Conference,  June  2005. 


122 


[41  j  Jones,  W.T.,  Nielsen,  E.J.  and  Park,  M.A.,  “Validation  of  3D  Adjoint  Based  Error 
Estimation  and  Mesh  Adaptation  for  Sonic  Boom  Prediction,”  Proceedings  of  the  44th 
AIAA  Aerospace  Sciences  Meeting  and  Exhibit,  Jan  2006. 

[42]  Bibb,  K.L.,  Gnoffo,  P.A.,  Park,  M.A.  and  Jones,  W.T.  ,  “Parallel,  Gradient-Based 
Anisotropic  Mesh  Adaptation  for  Re-entry  Vehicle  Configurations,”  Proceedings  of  the 
9th  AIAA/ASME  Joint  Therophysics  and  Heat  Transfer  Conference,  Jun  2006. 

[43]  Sankaran,  V.,  Sitaraman,  J.,  Wissink,  A.,  Datta,  A.,  Jayaraman,  B.,  Potsdam,  M., 
Mavriplis,  D.,  Yang,  Z.,  O’Brien,  D.,  Saberi,  H.,  Cheng,  R.,  Hariharan,  N.,  and  Strawn, 
R.,  “Application  of  the  Helios  Computational  Platform  to  Rotorcraft  Flowfields,”  Pro¬ 
ceedings  of  the  48th  AIAA  Aerospace  Sciences  Meeting  and  Exhibit,  Jan  2010. 

[44]  Wissink,  A.M.,  Kamkar.  S.,  Pulliam,  T.H.,  Sitaraman,  J.,  and  Sankaran,  V.,  “Cartesian 
Adaptive  Mesh  Refinement  for  Rotorcraft  Wake  Resolution,”  Proceedings  of  the  28th 
AIAA  Applied  Aerodynamics  Conference,  June  2010. 

[45]  Wissink,  A.,  Potsdam,  M.,  Sankaran,  V.,  Sitaraman,  J.,  Yang,  Z.,  and  Mavriplis,  D.  J., 
“A  Coupled  Unstructured-Adaptive  Cartesian  CFD  Approach  for  Hover  Prediction,” 
Proceedings  of  the  66th  Annual  Forum  of  the  American  Helicopter  Society,  May  2010. 

[46]  Stcrcnborg.  J.,  van  Zuijlcn,  A.,  and  Bijl,  H.,  “Solution  Based  Mesh  Adaptation  Applied 
to  Fluid  Structure  Interaction  Computations,”  47th  AIAA  Aerospace  Sciences  Meeting, 
Jan  2009. 

[47]  Alauzet,  F.  and  Olivier,  G.,  “Extension  of  Metric-Based  Anisotropic  Mesh  Adaptation 
to  Time-Dependent  Problems  Involving  Moving  Geometries,”  Proceedings  of  the  49th 
AIAA  Aerospace  Sciences  Meeting  and  Exhibit,  Jan  2011. 

[48]  Pirzadeh,  S.,  “Advanced  Unstructured  Grid  Generation  for  Complex  Aerodynamic  Ap¬ 
plications,”  AIAA  Journal ,  Vol.  48,  (5),  2010,  pp.  904-915. 

[49]  Noack.  R.,  Boger,  D.,  Kunz,  R.,  and  Carrica,  R,  “SUGGAR+- f:  An  Improved  General 
Overset  Grid  Assembly  Capability,”  Proceedings  of  the  19th  AIAA  Computational  Fluid 
Dynamics  Conference,  June  2009. 

[50]  Renaud,  T.,  O'Brien,  D.  M.,  Jr.,  Smith,  M.  J.,  and  Potsdam,  M.,  “Evaluation  of  Iso¬ 
lated  Fuselage  and  Rotor-Fuselage  Interaction  Using  CFD,”  Journal  of  the  American 
Helicopter  Society ,  Vol.  53,  (1),  2008,  pp.  3-17. 

[51]  Liggett,  N.  and  Smith,  M.  J.,  “Temporal  Convergence  Criteria  for  Time- Accurate 
Viscous  Simulations  of  Separated  Flows,”  Computers  &  Fluids ,  Vol.  66,  doi: 
10.1016/j.compfluid.2012.06.010,  2012,  pp.  140-156. 

[52]  Komerath,  N.,  “Deconstructing  Hub  Drag,  Part  I:  Experiments,”  Technical  Report 
ONR  Final  Report,  September  2013. 


123 


[53]  Anderson.  W.,  Rausch.  R.,  and  Bonhaus,  D.,  1  Implicit/Multigrid  Algorithms  for  Incom¬ 
pressible  Turbulent  Flows  on  Unstructured  Grids,”  Journal  of  Computational  Physics, 
Vol.  128,  (2),  1996. 

[54]  Chorin,  A.,  “A  Numerical  Method  for  Solving  Incompressible  Viscous  Flow  Problems,” 
Journal  of  Computational  Physics,  Vol.  135,  (2),  1997. 

[55]  Roe,  P.,  “Approximate  Riemann  Solvers,  Parameter  Vectors,  and  Difference  Schemes,” 
Journal  of  Computational  Physics,  Vol.  43,  (2),  1981. 

[56]  Menter,  F.,  “Two-Equation  Eddy-Viscosity  Turbulence  Models  for  Engineering  Appli¬ 
cations,”  AIAA  Journal,  Vol.  32,  (8),  1994. 

[57]  Sanchez-Rocha,  M.  and  Menon,  S.,  “The  Compressible  Hybrid  RANS/LES  Formula¬ 
tion  using  an  Additive  Operator,”  Journal  of  Computational  Physics,  Vol.  228,  2009, 
pp.  2037-2062. 

[58]  Sanchez-Rocha,  M.  and  Menon,  S.,  “An  order-of-magnitude  approximation  for  the  hy¬ 
brid  terms  in  the  compressible  hybrid  RANS/LES  governing  equations,”  Journal  of 
Turbulence,  Vol.  12,  2011,  pp.  1-22. 

[59]  Noack,  R.,  ‘DiRTlib:  A  Library  to  Add  an  Overset  Capability  to  Your  Flow  Solver,” 
Proceedings  of  the  17th  AIAA  Computational  Fluid  Dynamics  Conference,  June  2005. 

[60,  Renaud,  T.,  D.  M.  O'Brien,  .L,  Smith.  M.  J.,  and  Potsdam.  M.,  Evaluation  of  Iso¬ 
lated  Fuselage  and  Rotor-Fuselage  Interaction  Using  CFD,”  Journal  of  the  American 
Helicopter  Society,  Vol.  53,  (1),  2008,  pp.  3-17. 

[61]  Lynch,  C.  E.,  Advanced  CFD  Methods  for  Wind  Turbine  Analysis ,  Ph.D.  thesis,  Georgia 
Institute  of  Technology,  2011. 

[62]  Abras,  J.  N.,  Enhancement  of  Aeroelastic  Rotor  Airload  Prediction  Methods,  Ph.D. 
thesis,  Georgia  Institute  of  Technology,  2009. 

[63]  Lynch,  C.  E.  and  Smith,  M.  J.,  “Extension  and  Exploration  of  a  Hybrid  Turbulence 
Model  on  Unstructured  Grids,”  AIAA  Journal,  Vol.  49,  (11),  2011,  pp.  2585-2591. 
doi:  10.2514/  1.56296 

[64]  Liggett,  N.  and  Smith,  M.  J.,  “A  Study  of  the  Cap  Physics  of  Airfoils  with  Unsteady 
Flaps,”  AIAA  Journal  of  Aircraft,  Vol.  50,  (2),  doi:  10.2514/1.C032026,  2013,  pp.  643- 
650. 

[65]  Kravchenko,  A.  G.  and  Moin,  P.,  “Numerical  Studies  of  Flow  Over  a  Circular  Cylinder 
at  Reo  =  3900,”  Physics  of  Fluids,  Vol.  12,  (2),  2000,  pp.  403-417. 

[66]  Norberg,  C.,  Unpublished  data  that  was  digitized  from  Fig.  11  in  Ref.65. 


124 


[67]  Son,  J.  and  Hanratty,  T.  J.,  “Velocity  Gradients  at  the  Wall  for  Flow  Around  a  Cylinder 
at  Reynolds  Numbers  5  x  103  to  105,”  Journal  of  Fluid  Mechanics ,  Vol.  35,  (2),  1969, 
pp.  353-368. 

[68]  Ong,  L.  and  Wallace,  J.,  “The  Velocity  Field  of  the  Turbulent  Very  Near  Wake  of  a 
Circular  Cylinder,”  Experimental  Fluids ,  Vol.  20,  (6),  1996,  pp.  441-453. 

[69]  Mentcr,  F.,  “Two- Equation  Eddy- Viscosity  Turbulence  Models  for  Engineering  Appli¬ 
cations,”  AIAA  Journal,  Vol.  32,  (8),  1994,  pp.  1598-1605. 

[70]  Spalart,  P.  R.,  “Strategies  for  Turbulence  Modelling  and  Simulations,”  International 
Journal  of  Heat  and  Fluid  Flow,  Vol.  21,  (3),  2000,  pp.  252-263. 

[71 1  Menter,  F.,  “Review  of  the  Shear-Stress  Transport  Turbulence  Model  Experience  from 
an  Industrial  Perspective,"  International  Journal  of  Computational  Fluid  Dynamics , 
Vol.  23,  (4),  2009,  pp.  305-316. 

[72]  Smirnov,  P.  E.  and  Menter,  F.  R.,  “Sensitization  of  the  SST  Turbulence  Model  to 
Rotation  and  Curvature  by  Applying  the  Spalart-Shur  Correction  Term,”  Journal  of 
Turbomachinery,  Vol.  131,  (4),  2009,  pp.  1-8. 

[73]  “tecplot  360  2013  User’s  Manual.  Release  1,”  Technical  report,  Teeplot,  Inc.,  Bellevue, 
WA,  2013. 

[74]  “Ficldvicw  13  User’s  Manual,”  Technical  report,  Intelligent  Light,  2013. 

[75]  “MATLAB  R2013b  User’s  Manual,”  Technical  report,  MathWorks,  Inc.,  Natick,  Mass- 
chusetts,  2013. 

[76 1  Bridgcman,  J.  and  Lancaster,  G.  “Physics-Based  Analysis  Methodology  for  Hub  Drag 
Prediction,”  Proceedings  of  the  66th  Annual  Forum  of  the  American  Helicopter  Society, 
May  2010. 

[77]  Shcnoy,  R.,  Smith,  M.  J.,  and  Park,  M.  A.,  “Unstructured  Overset  Mesh  Adaptation 
with  Turbulence  Modeling  for  Unsteady  Aerodynamic  Interactions,”  AIAA  Journal  of 
Aircraft,  in  press. 

78]  Park,  M.  A.,  Anisotropic  Output-Based  Adaptation  with  Tetrahedral  Cut  Cells  for  Com¬ 
pressible  Flows,  Ph.D.  thesis,  Massachusetts  Institute  of  Technology,  2008. 

[79]  Kamkar,  S.  J.,  Mesh  Adaption  Strategies  for  Vortex- Dominated  Flows,  Ph.D.  thesis, 
Stanford  University,  2011 

[80]  Venditti,  D.  A.,  Grid  Adaptation  for  Functional  Outputs  of  Compressible  Flow  Simula¬ 
tions,  Ph.D.  thesis,  Massachusetts  Institute  of  Technology,  2002. 

[81  Brand,  A.  G.,  An  Experimental  Investigation  of  the  Interaction  Between  a  Model  Rotor 
and  Airframe  in  Forward  Flight,  Ph.D.  thesis,  Georgia  Institute  of  Technology,  1989. 


125 


[82]  Mavris,  D.  N.,  An  Analytical  Method  for  the  Prediction  of  Unsteady  Rotor /Airframe 
Interactions  in  Forward  Flight,  Ph.D.  thesis,  Georgia  Institute  of  Technology,  1988, 
revised  1991. 

[83]  O’Brien,  D.,  and  Smith,  M.,  “Analysis  of  Rotor-Fuselage  Interactions  Using  Various 
Rotor  Models,”  Proceedings  of  the  43rd  AIAA  Aerospace  Sciences  Meeting  and  Exhibit, 
Jan  2005. 

[84]  Lee,  J.  and  Kwon,  O.J.,  “Predicting  Aerodynamic  Fuselage  Interactions  by  Using  Un¬ 
structured  Meshes,”  Transactions  of  the  Japanese  Society  for  Aeronautical  and  Space 
Sciences,  Vol  44,  (146),  2002 

[85]  O’Brien,  D.  M.,  Analysis  of  Computational  Modeling  Techniques  for  Complete  Rotor- 
craft  Configurations ,  Ph.D.  thesis,  Georgia  Institute  of  Technology,  Atlanta,  Georgia, 
2006. 

[86]  Park,  M.A.  and  Carlson.  J.R.,  “Turbulent  Output-Based  Anisotropic  Adaptation,”  Pro¬ 
ceedings  of  the  48th  AIAA  Aerospace  Sciences  Meeting  and  Exhibit,  Jan  2010. 

[87]  Mineck,  R.  and  Althoff  Gorton,  S.,  “Steady  and  Periodic  Pressure  Measurements  on 
a  Generic  Helicopter  Fuselage  Model  in  the  Presence  of  a  Rotor,"  Technical  Report 
NASA/TM-2000-210286,  June  2000. 

[88]  Ghee,  T.A.  and  Elliott,  J.W.,  “The  Wake  of  a  Small-Scale  Rotor  in  Forward  Flight 
Using  Flow7  Visualization,”  Journal  of  the  American  Helicopter  Society ,  Vol.  40,  (3), 
1995. 

[89]  Park,  Y.M.,  Nam,  H.J,  and  Kwon,  O..J.,  “Simulation  of  Unsteady  Rotor-Fuselage  In¬ 
teractions  Using  Unstructured  Adaptive  Meshes,”  Journal  of  the  American  Helicopter 
Society,  Vol.  51,  (2),  2006. 

[90]  Kenyon,  A.R.  and  Brown,  R.E.,  "Wake  Dynamics  and  Rotor-Fuselage  Aerodynamic 
Interactions,”  Journal  of  the  American  Helicopter  Society ,  Vol.  54,  (1),  2009. 

[91]  Raghav,  V.,  Shenoy.  R.,  Smith,  M.  J.,  and  Komerath,  N.  M.,  “Investigation  of  Drag 
and  Wake  Turbulence  of  a  Rotor  Hub,”  Aerospace  Science  and  Technology ,  Vol.  28,  (1), 
doi:  10.1016/j.ast.2012.10.012,  2013,  pp.  164-175. 

[92]  Lynch,  C.  E.  and  Smith.  M.  J.,  “Hybrid  RANS-LES  Turbulence  Models  on  Unstructred 
Grids,”  Proceedings  of  the  38th  Fluid  Dynamics  Conference  and  Exhibit,  June  2008. 

[93]  Shenoy,  R.,  Holmes,  M.,  Smith,  M.  J.,  and  Komerath,  N.,  “Scaling  Evaluations  on  the 
Drag  of  a  Hub  System,”  Journal  of  the  American  Helicopter  Society,  Vol.  58,  (3),  July 
2013,  pp.  1-13. 

[94]  Schlichting,  II..  Boundary  Layer  Theory,  7th  edition,  McGraw-Hill,  New7  York,  NY, 
1979. 


126 


[95]  Hoerner,  S.  F.,  Fluid  Dynamic  Drag ,  Sighar  Hoerner,  1965. 

[96]  Shenoy,  R.,  Smith,  M.  J.,  and  Komerath,  N.,  “Computationally  Deconstructing  Drag 
From  Model  To  Full  Scale,”  Technical  Report  Vertical  Lift  Consortium  Final  Report, 
WBS  201 1-B-l  1-07. 3-A1,  February  2013. 

[97]  Zdravkovich,  M.  M.,  Flow  Around  Circular  Cylinders,  Vol.  2:  Applications,  Oxford 
Science  Publications,  New  York,  NY,  2003. 

[98]  Sakamoto,  H.  and  Arie,  M.,  “Vortex  shedding  from  a  rectangular  prism  and  a  circular 
cylinder  placed  vertically  in  a  turbulent  boundary  layer,”  Journal  of  Fluid  Mechanics, 
Vol.  126,  (-1),  1983,  pp.  147-165. 

[99]  Bousman,  W.  G.  and  Kufeld,  R.  M.,  “UH-60A  Airloads  Catalog,”  Technical  Report 
NASA  TM-2005-212827 / AFDD  TR-05-003,  August  2005. 


127 


REPORT  DOCUMENTATION  PAGE 

Form  Approved 

OMB  No.  0704-0188 

Public  reporting  burden  (or  this  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  data  sources, 

gathering  and  maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  information  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection 

of  information,  including  suggestions  for  reducing  this  burden  to  Washington  Headquarters  Service,  Directorate  for  Information  Operations  and  Reports, 

1215  Jefferson  Davis  Highway  Suite  1204,  Arlington  VA  22202-4302,  and  to  the  Office  of  Management  and  Budget, 

Paperwork  Reduction  Project  (0704-0188)  Washington,  DC  20503 

PLEASE  DO  NOT  RETURN  YOUR  FORM  TO  THE  ABOVE  ADDRESS. 

1.  REPORT  DATE  (DD-MM-YYYY) 

30-09-2013 

2.  REPORT  TYPE 

Final  Technical  Report  (Part  2  of  2) 

3.  DATES  COVERED  (From  -  To) 

01-07-2009  to  30-06-2013 

4.  TITLE  AND  SUBTITLE 

Deconstructing  Hub  Drag 

Part  II:  Computational  Development  and  Analysis 

5a.  CONTRACT  NUMBER 

N0001409-1-1 01 9 


5c.  PROGRAM  ELEMENT  NUMBER 


6.  AUTHOR(S) 

Rajiv  Shenoy  and  Marilyn  J.  Smith 


5d.  PROJECT  NUMBER 


5e.  TASK  NUMBER 


5f.  WORK  UNIT  NUMBER 


7.  PERFORMING  ORGANIZATION  NAME(S)  AND  ADDRESS(ES) 
Georgia  Tech  Research  Corporation 
Daniel  Guggenheim  School  of  Aerospace  Engineering 
Georgia  Institute  of  Technology 
Atlanta,  GA  30332-0150 


8.  PERFORMING  ORGANIZATION 
REPORT  NUMBER 


9.  SPONSORING/MONITORING  AGENCY  NAME(S)  AND  ADDRESS(ES) 
Office  of  Naval  Research 
875  North  Randolph  Street 
Arlington,  VA  22203-1995 


10.  SPONSOR/MONITOR'S  ACRONYM(S) 


11.  SPONSORING/MONITORING 
AGENCY  REPORT  NUMBER 


12.  DISTRIBUTION  AVAILABILITY  STATEMENT 

Distribution  approved  for  public  release;  distribution  is  unliminited 


13.  SUPPLEMENTARY  NOTES 

Part  I  of  this  report  on  the  experimental  research  reported  separately  by  co-PI  Dr.  N.  Komerath 


14.  ABSTRACT 

Development,  validation  and  demonstration  of  a  computational  approach  to  model  complex  rotorcraft  hubs,  their 
components  and  their  wakes  for  both  static  and  rotating  operation  has  been  developed.  A  priori  computations  were 
compared  with  experimental  results  to  yield  high  confidence  in  the  computational  analysis.  The  hub  was  analyzed  for  its 
total  and  component  drag,  including  interference  drag.  Components  include  hub  plates,  shanks,  pitch  links,  scissors, 
and  hardware.  Wake  turbulent  spectra  were  obtained  for  the  near  and  far  field  wakes,  and  wake  velocities  in  the  near 
wake.  Scaling,  influence  of  the  fuselage,  and  physics  of  the  wake  were  all  explored.  With  the  highly  accurate 
computational  approach  developed,  design  and  analysis  of  hubs  can  be  conducted  using  computations  directly  on  full- 
scale  configurations. 


15.  SUBJECT  TERMS 

Hub  Drag,  Computational  Fluid  Dynamics,  CFD,  Large  Eddy  Simulations,  LES,  turbulence  modeling,  grid  adaptation, 
rotational  effects,  turbulent  wake  spectra,  drag,  interference  drag, 


16.  SECURITY  CLASSIFICATION  OF: 

17.  LIMITATION  OF 
ABSTRACT 

Unlimited 

18.  NUMBER 

OF  PAGES 

127 

19a.  NAME  OF  RESPONSIBLE  PERSON 

Prof.  Marilyn  J.  Smith 

a.  REPORT 

Unclassified 

b.  ABSTRACT 

Unclassified 

c.  THIS  PAGE 

Unclassified 

19b.  TELEPONE  NUMBER  ( Include  area  code) 

404-894-3065 

Standard  Form  298  (Rev.  8-98) 
Prescribed  by  ANSI-Std  Z39-18 


